{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# VinDsl.jl: Fast and furious statistical modeling\n", " \n", "\n", "\n", "\n", "
\n", "John Pearson \n", "P[λ]ab \n", "Duke Institute for Brain Sciences\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Following along\n", "VinDsl currently makes use of some features of Distributions.jl that are not yet available on master, as well as the latest release of PDMats.jl. You will need to checkout the `jp/autodiff` branch of Distributions:\n", "```julia\n", "Pkg.clone(\"https://github.com/jmxpearson/VinDsl.jl\")\n", "Pkg.update() # if you don't have latest Distributions\n", "Pkg.checkout(\"Distributions\", \"jp/autodiff\")\n", "```" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true, "slideshow": { "slide_type": "slide" } }, "source": [ "# About me\n", "\n", "\n", "\n", "- computational neuroscience lab at Duke\n", "- using Julia for about a year\n", "- member of JuliaStats organization" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Our problem \n", " \n", "
\n", "\n", "
By AKS.9955, via Wikimedia Commons
\n", "
\n", "\n", "
\n", "\n", "
https://praneethnamburi.wordpress.com/2015/02/05/simulating-neural-spike-trains/
\n", "
\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Our problem\n", "\n", "- model responses to known features (like GLM)\n", "- infer latent features\n", "- do this scalably for large populations of neurons\n", "- use (generative) Bayesian models that account for uncertainty" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "# One solution: Sampling\n", "\n", "- Markov Chain Monte Carlo (MCMC) methods come with guarantees about correctness\n", "- lots of packages (Lora.jl, Mamba.jl, Stan.jl)\n", "- **But** sampling does not scale well to millions of observations and parameters" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Another solution: The Max Power Way\n", "\n", "" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Variational Bayesian Inference (VB)\n", "\n", "Generative model for data: $p(y|\\theta)p(\\theta)$ \n", "Actual posterior: $p(\\theta|y) = p(y|\\theta)p(\\theta)/p(y)$ \n", "Approximate posterior: $q(\\theta)$ \n", "\n", "Maximize **E**vidence **L**ower **Bo**und (ELBO) wrt $\\theta$:\n", "\n", "$$\n", "\\mathcal{L} = \\log p(y) \\ge -KL\\left(q \\middle\\| p\\right) = \\mathbb{E}_q[\\log p(y|\\theta)] + \\mathcal{H}[q]\n", "$$\n", "\n", "- $\\mathbb{E}_q[\\log p(y|\\theta)]$ measures goodness of fit (data likely under approximate posterior)\n", "- $\\mathcal{H}[q]$ is the entropy (favors less certain models)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Why VB?\n", "\n", "- Scales well\n", "- Can use well-studied optimization techniques\n", "\n", "# Drawbacks:\n", "- !@$*&# hard to code\n", "- Can't quickly spec out a model like with Stan or JAGS/BUGS\n", "- Traditionally, requires that distributions be conjugate, requires doing lots of algebra" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# But VB is exploding!\n", "\n", "- stochastic variational inference (SVI): [Hoffman et al.](http://dl.acm.org/citation.cfm?id=2502622)\n", "- black box variational inference (BBVI): [Ranganath et al.](http://arxiv.org/abs/1401.0118)\n", "- control variates: [Paisley et al.](http://arxiv.org/abs/1206.6430)\n", "- local expectation gradients (LEG): [Titsias and Lázaro-Gredilla](http://papers.nips.cc/paper/5678-local-expectation-gradients-for-black-box-variational-inference)\n", "- neural variational inference (NVIL): [Mnih and Gregor](http://arxiv.org/abs/1402.0030)\n", "- variational autoencoders [Kingma and Welling](http://arxiv.org/abs/1312.6114), [Rezende et al.](http://arxiv.org/abs/1401.4082) " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# What's out there\n", "\n", "- research code: individual algorithms, proof of concept \n", "- [Stan](http://mc-stan.org/) \n", " - **Pros**: Robust, actively developed, good with stats, accessible from Julia\n", " - **Cons**: variational methods still experimental, C++\n", " " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "- [Edward](https://github.com/blei-lab/edward)\n", " - **Pros**: Python, multiple backends, from Stan and VB core developers\n", " - **Cons**: Python, very new\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# What's out there\n", "- [Theano](http://deeplearning.net/software/theano/) & [TensorFlow](https://www.tensorflow.org/)\n", " - **Pros**: well-tested, stable, well-engineered backends\n", " - **Cons**: complex, C++, not very stats-aware\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# What do we want?\n", "- write math, get code — a domain-specific language (DSL)\n", "- easily generalize to different numbers of indices, structures\n", "- only weakly opinionated about model structure or inference\n", "- model code should be *hackable*\n", " - easy to use prefab pieces\n", " - not hard to write custom VB tricks\n", " - fast prototyping\n", "- no (or minimal) algebra\n", "- automatic gradients" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true, "slideshow": { "slide_type": "slide" } }, "source": [ "Introducing...\n", "![](http://www.joblo.com/newsimages1/vin.diesel_1920x1200_961)\n", "## VinDsl.jl: Fast and furious variational inference" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# What's the goal?\n", "- quick prototyping\n", "- targeted at researchers\n", "- \"internal\" DSL\n", "- loosely-coupled parts\n", "- \"consenting adults\" philosophy" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "# Today: two prototypes\n", "- \"Classical\" models (conjugate + some optimization)\n", "- ADVI \"Black Box\" models" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Model 1\n", "\n", "### Finally, some code!" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false, "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "WARNING: Method definition logpdf(Distributions.Gamma, Real) in module Distributions at /Users/jmxp/.julia/v0.4/Distributions/src/univariates.jl:321 overwritten in module VinDsl at /Users/jmxp/.julia/v0.4/VinDsl/src/distributions/expfam.jl:42.\n", "WARNING: Method definition logpdf(Distributions.Poisson, Int64) in module Distributions at /Users/jmxp/.julia/v0.4/Distributions/src/univariates.jl:321 overwritten in module VinDsl at /Users/jmxp/.julia/v0.4/VinDsl/src/distributions/expfam.jl:64.\n" ] } ], "source": [ "using Distributions\n", "using VinDsl" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Model structure:\n", "### Main idea: Factor graphs\n", "- idea from Dahua Lin in [this talk](http://people.csail.mit.edu/dhlin/jubayes/julia_bayes_inference.pdf)\n", "\n", "- Nodes: arrays of distributions\n", "- Factors $\\leftrightarrow$ terms in variational objective\n", " - but not locked in to graphical model structure!" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "dims = (20, 6)\n", "\n", "μ[j] ~ Normal(zeros(dims[2]), ones(dims[2]))\n", "τ[j] ~ Gamma(1.1 * ones(dims[2]), ones(dims[2]))\n", "μ0[j] ~ Const(zeros(dims[2]))\n", "\n", "y[i, j] ~ Const(rand(dims));" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "notes" } }, "source": [ "Nodes: under the hood\n", "nodes define the q/approximate posterior/recognition model\n", "~ defines a node\n", "can use any distribution defined in the Distributions package\n", "code parses the left and right-hand sides\n", "indices on left get tracked and assigned to dimensions of parameter arrays\n", "code is rewritten as a call to a node constructor" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "f = @factor LogNormalFactor y μ τ;" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "or..." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false, "slideshow": { "slide_type": "-" } }, "outputs": [ { "data": { "text/plain": [ "1-element Array{VinDsl.Factor{N},1}:\n", " VinDsl.LogNormalFactor{2}\n", "x: VinDsl.ConstantNode{Float64} y[i,j]\n", "μ: VinDsl.RandomNode{Normal} μ[j]\n", "τ: VinDsl.RandomNode{Gamma} τ[j]\n" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "@pmodel begin\n", " y ~ Normal(μ, τ)\n", "end" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "New factor types can be defined with yet another macro:" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "```julia\n", "@deffactor LogNormalFactor [x, μ, τ] begin\n", " -(1/2) * ((E(τ) * ( V(x) + V(μ) + (E(x) - E(μ))^2 ) + log(2π) + Elog(τ)))\n", "end\n", "\n", "@deffactor LogGammaFactor [x, α, β] begin\n", " (E(α) - 1) * Elog(x) - E(β) * E(x) + E(α) * E(β) - Eloggamma(α)\n", "end\n", "```" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "- Uses a \"mini-language\" with `E(x)` $\\equiv \\mathbb{E}[X]$, `V(x)` $\\equiv \\textrm{cov}[X]$, etc. \n", "- Again, no need to track indices\n", " - multivariate distributions (Dirichlet, MvNormal) are automatically multivariate in these expressions\n", "- `VinDsl` generates a `value(f)` function that handles indices appropriately and sums over the dimensions of the array" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Models are just factor graphs\n", "Let's do a simple conjugate model: \n", "p-model:\n", "$$\n", "\\begin{align}\n", "\\mu_j &\\sim \\mathcal{N}(0, 1) \\\\\n", "\\tau_j &\\sim \\mathrm{Gamma}(1.1, 1) \\\\\n", "y_{ij} &\\sim \\mathcal{N}(\\mu_j, \\tau_j)\n", "\\end{align}\n", "$$\n", "\n", "q-model:\n", "$$\n", "\\begin{align}\n", "\\mu &\\sim \\mathcal{N}(m, t) \\\\\n", "\\tau &\\sim \\mathrm{Gamma}(\\alpha, \\beta) \\\\\n", "\\end{align}\n", "$$" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": true, "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "dims = (20, 6)\n", "\n", "# note: it won't matter how we initialize here\n", "μ[j] ~ Normal(zeros(dims[2]), ones(dims[2]))\n", "τ[j] ~ Gamma(1.1 * ones(dims[2]), ones(dims[2]))\n", "μ0[j] ~ Const(zeros(dims[2]))\n", "τ0[j] ~ Const(2 * ones(dims[2]))\n", "a0[j] ~ Const(1.1 * ones(dims[2]))\n", "b0[j] ~ Const(ones(dims[2]))\n", "\n", "y[i, j] ~ Const(rand(dims))\n", "\n", "@pmodel begin\n", " y ~ Normal(μ, τ)\n", " μ ~ Normal(μ0, τ0)\n", " τ ~ Gamma(a0, b0)\n", "end\n", "\n", "m = VBModel([μ, τ, μ0, τ0, a0, b0, y], pmodel_factors);" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true, "slideshow": { "slide_type": "slide" } }, "source": [ "# What's going on here?" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "Dict{VinDsl.Node,Symbol} with 7 entries:\n", " VinDsl.RandomNode{Gamma… => :conjugate\n", " VinDsl.ConstantNode{Flo… => :constant\n", " VinDsl.ConstantNode{Flo… => :constant\n", " VinDsl.ConstantNode{Flo… => :constant\n", " VinDsl.ConstantNode{Flo… => :constant\n", " VinDsl.RandomNode{Norma… => :conjugate\n", " VinDsl.ConstantNode{Flo… => :constant" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "m.update_strategy" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [], "source": [ "update!(m)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Conjugacy\n", "- right now `VinDsl` goes out of its way to handle conjugacy between nodes\n", "- conjugate relationships not automatically detected, but easy to define\n", "- `@defnaturals` returns expected sufficient statistics from a factor for a given target distribution" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "```julia\n", "@defnaturals LogNormalFactor μ Normal begin\n", " Ex, Eτ = E(x), E(τ)\n", " (Ex * Eτ, -Eτ/2)\n", "end\n", "```" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Index Bookkeeping\n", "- nodes have associated indices\n", "- factors know which indices go with which nodes, which indices to sum over\n", " - inner indices belong to, e.g., elements of a multivariate normal (should not be separated)\n", " - outer indices correspond to replicates of \"atomic\" variables " ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "So this is easy: `i` is inner:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [], "source": [ "d = 5\n", "μ[i] ~ MvNormalCanon(zeros(d), diagm(ones(d)))\n", "Λ[i, i] ~ Wishart(float(d), diagm(ones(d)));" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "But here, `i` is inner for $\\mu$ but not for $\\tau$. In any factor combining these two, $\\tau$ will be treated like a vector because it matches an inner index for some node:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": true }, "outputs": [], "source": [ "μ[i] ~ MvNormalCanon(zeros(d), diagm(ones(d)))\n", "τ[i] ~ Gamma(1.1 * ones(d), ones(d));" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Model 2: ADVI\n", "[(Automatic Differentiation Variational Inference)](http://arxiv.org/abs/1603.00788)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true, "slideshow": { "slide_type": "slide" } }, "source": [ "# Two major ideas\n", "1. We just need to define an ELBO (or *approxmate* an ELBO)\n", "1. Any unimodal distribution is *approximately* a normal with constrained support" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "# Implementation\n", "- Approximate ELBO by sampling from normal, transforming to constrained variables\n", "- Let automatic differentiation handle the gradient calculation\n", "- Do gradient ascent" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Hey, we have autodiff!\n", "![http://www.juliadiff.org/](juliadiff.png)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# But there's always a problem:\n", "- ForwardDiff defines `Dual <: Real`\n", "- But Distributions.jl doesn't allow, e.g., `Dual` $\\mu$ and $\\sigma$ for `Normal`\n", "- Then a lot of work..." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Today:\n", "- Some distributions allow arbitrary parameter types:\n", " - using PDMats:\n", " - MvNormal, MvNormalCanon, GenericMvTDist\n", " - Wishart, InverseWishart\n", " - other:\n", " - Dirichlet, Normal, NormalCanon, Gamma, InverseGamma, Poisson\n", "- for examples see [here](https://github.com/jmxpearson/distribution_diff_tests)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# This month:\n", "- ~40 more distributions — almost all univariates\n", "- PR in progress from @Halmoni100 in my lab\n", "- some special functions in StatsFuns and Base still assume Float64" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Defining a simple model\n", "$$\\begin{align}\n", "\\sigma_\\eta &\\sim \\mathrm{Gamma}(a_\\eta, b_\\eta) \\\\\n", "a_u &\\sim \\mathcal{N}(\\mu_a, \\sigma_a) \\\\\n", "\\eta_{tu} &\\sim \\mathcal{N}(a_u, \\sigma_\\eta) \\\\\n", "N_{tu} &\\sim \\mathrm{Poisson}(e^{\\eta_{tu}}) \\\\\n", "\\end{align}$$" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "text/plain": [ "ELBO (generic function with 1 method)" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "@ELBO begin\n", " @advi_declarations begin\n", " a::Real()[U]\n", " σ::Positive()\n", " η::Real()[T, U] \n", " end\n", " \n", " @advi_model begin\n", " for u in 1:U\n", " a[u] ~ Normal(log(15.), 0.1)\n", " end\n", " \n", " σ ~ Gamma(1, 1)\n", "\n", " for t in 1:T\n", " for u in 1:U\n", " η[t, u] ~ Normal(a[u], σ)\n", " spikes[t, u] ~ Poisson(exp(η[t, u]))\n", " end\n", " end\n", " end\n", "end" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false, "slideshow": { "slide_type": "skip" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/Users/jmxp/anaconda/lib/python2.7/site-packages/matplotlib/font_manager.py:273: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.\n", " warnings.warn('Matplotlib is building the font cache using fc-list. This may take a moment.')\n" ] } ], "source": [ "using PyPlot\n", "T = 100 # number of time steps\n", "U = 10 # units\n", "\n", "baseline = 10.\n", "baseline_sd = 0.10\n", "log_bl = log(baseline)\n", "unit_bl = log_bl + baseline_sd * randn(U);\n", "\n", "fr_log = unit_bl' .+ zeros(T, U) # firing rate\n", "\n", "σ_η = 0.05\n", "\n", "eta = fr_log + σ_η * randn(T, U)\n", "fr = exp(eta);\n", "\n", "spikes = Array{Int}(size(fr)...)\n", "for i in eachindex(fr)\n", " spikes[i] = rand(Poisson(fr[i]))\n", "end" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Add some data..." ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABPwAAADGCAYAAABGk54uAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X1wVuWd//FPIs8mOAlIE9kCZhN8ADNJ6DJsaQbcOqOhEp1ZEnxodQnFXZ0QarvWaVzUlULr1K0udAYKU0Jb7Cqx4JZKSovVMSrCAAFnKkrSHXkICSyQCpkQCHD//tgftCkE87nIfWMO79cM/5z7e3Jd51yP58v9kBSLxWICAAAAAAAAEAnJl7sCAAAAAAAAAHoOCT8AAAAAAAAgQkj4AQAAAAAAABFCwg8AAAAAAACIEBJ+AAAAAAAAQISQ8AMAAAAAAAAihIQfAAAAAAAAECEk/AAAAAAAAIAIIeEHAAAAAAAARAgJPwAAAAAAACBCLkvC7+TJk3r88cc1fPhwDRo0SBMmTNCGDRsuR1UABNqyZYvKy8s1duxYpaSkaOTIkZo+fbrq6+vPi/3www91xx13KDU1VUOGDNEDDzygQ4cOXYZaAwg1f/58JScnKzc397zXGONA77Nt2zYVFxdryJAhuvrqq3XLLbfoRz/6UacYxjbQ+zQ0NOiee+7R5z//eV199dW66aabNG/ePB0/frxTHOMbiL6kWCwWS3Sh9957r1avXq1HH31U2dnZWrFihTZv3qw333xTX/ziFxNdHQABSkpK9O6776qkpES5ublqbm7WokWL1Nraqk2bNunmm2+WJDU2NiovL09paWmaM2eOjh07ph/84AcaOXKkNm/erD59+lzmKwHwaRobG3XDDTcoOTlZo0aN0vvvv9/pNcY40Lv89re/VXFxsQoKCjR9+nSlpKToj3/8o86cOaPvf//7khjbQG+0b98+3XLLLUpLS9O//Mu/KD09XRs3blRVVZXuuusurVmzRhLjG7hixBJs06ZNsaSkpNgPf/jDc8fa29tj2dnZsYkTJya6OgACbdy4MdbR0dHpWH19fWzAgAGxr33ta+eOPfzww7Grr746tm/fvnPHNmzYEEtKSootW7YsYfUFEG769Omx2267LTZ58uTYLbfc0uk1xjjQuxw9ejSWkZERmzZt2kXjGNtA7zN//vxYcnJybOfOnZ2OP/jgg7Hk5OTYn/70p1gsxvgGrhQJ/0jvK6+8oj59+mjWrFnnjvXv318zZ87Uxo0b1djYmOgqAQgwYcKE8/73Lzs7W2PGjNHOnTvPHVu9erXuvPNODR8+/NyxL3/5yxo9erRWrVqVsPoCCPPWW29p9erVeuGFFy74OmMc6F1efPFFHTx4UPPnz5cktbW1KXaBD/wwtoHe59ixY5KkYcOGdTqekZGh5ORk9evXTxLjG7hSJDzht337do0ePVopKSmdjo8fP/7c6wB6rwMHDmjo0KGSpP379+vgwYP6whe+cF7c+PHjVVdXl+jqATCcOXNGFRUVmjVrlsaMGXPe64xxoPd5/fXXNXjwYO3du1c33nijUlJSNHjwYD3yyCM6ceKEJMY20FtNnjxZsVhMZWVl2rFjh/bt26eXX35ZS5Ys0Zw5czRw4EDGN3AFSXjCr6mpSZmZmecdz8zMVCwW0/79+xNdJQA9ZOXKlWpsbNQ999wj6f/Gu6Qux/yRI0fU0dGR0DoC6L7Fixdrz549mjdv3gVfZ4wDvU99fb06Ojp01113qaioSKtXr9bMmTO1ZMkSlZWVSWJsA73V7bffrnnz5ul3v/ud8vPzNWLECN13332qqKjQc889J4nxDVxJEv5tnMePH1f//v3POz5gwIBzrwPofT788EOVl5dr4sSJeuCBByT9eTx/2pjv27dv4ioKoFuOHDmip556Sk8++aTS09MvGMMYB3qf1tZWHT9+XA8//LCef/55SdLdd9+tEydOaOnSpXrmmWcY20AvNmrUKE2aNEnTpk1Tenq6XnvtNc2fP18ZGRl65JFHGN/AFSThCb+BAwee+7jAX2pvbz/3OoDe5cCBA/rKV76itLQ0VVdXKykpSdKfxzNjHuh9nnjiCQ0ZMkTl5eVdxjDGgd7n7Jg8+278s+677z79+Mc/1saNG3XTTTdJYmwDvc1LL72khx56SA0NDefewXf33Xfr9OnTevzxx3XvvfeydgNXkIR/pDczM/Pc24j/0tlj1113XaKrBOASHD16VHfccYeOHj2q3/zmN8rIyDj32tmNRldjPj09nf89BD6DGhoatGzZMlVUVKixsVG7d+/Wxx9/rPb2dnV0dGj37t1qaWlhjAO90Nm99uc+97lOx89+yT9jG+i9Fi9erIKCgvM+rltcXKy2tjbV1dUxvoErSMITfnl5edq1a5daW1s7HX/vvfeUlJSkvLy8RFcJQKATJ07ozjvvVENDg1577TXdcMMNnV6/7rrrdO2112rLli3nnbt582bGO/AZ1djYqFgspoqKCl1//fW6/vrrlZWVpU2bNumjjz5SVlaW5s2bxxgHeqFx48ZJ+r9x/pfOfo/2sGHDGNtAL3XgwAGdPn36vOMdHR2KxWI6deoU4xu4giQ84Tdt2jSdOnVKS5cuPXfs5MmTWrFihSZMmNDpp8EBfHadOXNGpaWl2rRpk1555ZVzv7T91/7xH/9Rv/71rzs9WLz++uvatWuXSktLE1VdAIaxY8dqzZo1WrNmjV599dVz/8aMGaORI0fq1Vdf1cyZMyUxxoHeprS0VLFYTD/5yU86HV+2bJn69u2rSZMmSWJsA73R6NGjVVdXp4aGhk7Hf/GLX+iqq65Sbm6uJMY3cKVIisVisUQXOn36dL366qv6xje+oezsbK1YsUJbtmzR73//e02cODHR1QEQ4Bvf+IYWLlyo4uJilZSUnPf6/fffL0nat2+fCgoKdM0112jOnDk6duyYnnvuOY0YMUKbN2/mIwNAL3Lrrbfq8OHDev/9988dY4wDvc/Xv/51VVVVqaSkRJMmTdIbb7yhX/7yl6qsrDz3q9yMbaD3qa2t1Ze//GWlp6ervLxcQ4YM0dq1a7V+/XrNmjVLS5YskcT4Bq4UlyXhd/LkSc2dO1crV65US0uLcnNz9d3vfle33XZboqsCINCtt96qt956q8vX//LjBDt37tQ3v/lNvf322+rXr5/uvPNOPffcc7r22msTUVUAPeTWW2/VkSNHtGPHjk7HGeNA73L69GktWLBAVVVV2r9/v0aOHKny8nLNnj27UxxjG+h9tmzZoqefflp1dXU6fPiwrr/+ev3TP/2THnvsMSUn//kDfoxvIPouS8IPAAAAAAAAQHwk/Dv8AAAAAAAAAMQPCT8AAAAAAAAgQkj4AQAAAAAAABFCwg8AAAAAAACIEBJ+AAAAAAAAQISQ8AMAAAAAAAAihIQfAAAAAAAAECF9ElHIoUOHtH79eo0aNUoDBw5MRJEAAAAAAABAZBw/flwff/yxbr/9dg0dOvSisXbC7+TJk5o7d65WrlyplpYW5ebm6rvf/a5uu+22Ls9Zv369vvrVr7pFAQAAAAAAAPgLK1eu1P3333/RGDvh9+CDD2r16tV69NFHlZ2drRUrVmjKlCl688039cUvfvGC54waNUqSNGLECA0YMKDTa42NjRo+fPh559xwww1u1VRQUGDF//d//7ddhluve++914p/+umnrXhJmjhxohX/+c9/3i6jrq7Oin/nnXes+Pz8fCte0qdms//aoEGD7DLccxJxb934kOt2r6OrPrh06VI99NBDF3zt8OHDVhltbW1W/KFDh6x4ye+HH330kV2GW68RI0ZY8W7/kPzrDunn7pzg9g/Jvw53bHz44YdWvOS3X0gZ7nW4deqq7Xbu3Kmbbrrpgq+5feTGG2+04iVpz5499jkOd84JOcdtixDuPBWyHrv91t2vSdKPfvQjK764uNguw20Pt71Dxre7xwtZAy40h2zYsOGi/5nvCFmPZ86cacWXl5fbZbjt7fapkPH97W9/24p375Mk/e53v7PiQ+YEd35217GQ+fmuu+6y4l966SW7jJD9kSNkLF1oDejq+VsKa+94CxlL7nx7/Phxuwy3vRPx7OqOb3eNkfx7G7LHc6/dfc4ImUPc8ReSz+runq29vV179uw5l2e7GCvht3nzZr388sv6j//4Dz366KOSpK997WsaO3asvv3tb+vtt9++4HlnP8Y7YMCA8xrvqquuumCDugNCkkaOHGnFp6Sk2GUMGzbMih8zZowVH1KnribsrmRlZdllNDY2WvF/ndj9NCHtnZGRYcWnpqbaZbjtkZ2dbZfR1NRkxbt1Crlut593dd1XX311l6+519Ha2mrF9+njf2OB235unSQpOdn76lR3M/Hxxx9b8SFlhPTz+vp6K/7MmTN2Ge51uH3wyJEjVrzk1ylkM+9eh1una6655oLH+/bt2+Vr7hzirt/S/33iIJ5Cxrd7jrt+hzh48KAVH/LA6vbbkPbu27evFe/uESS/Xm57h4xvd8/m7tekC88h/fv3D7qHPcXdP7v7Tqnrua0rblvk5ORY8ZJ/HSFjyb3ukHnKfZiO955Q8h++09LS7DIu55jpyt69e8871tXzt+Sv34mQiHUppE+57Z2ZmWnFh+QF3PGdiHsbUoZ77R988IEV7+4pJH9eC8lvXGi8Xkx3vi7PevJ85ZVX1KdPH82aNevcsf79+2vmzJnauHFj0CYDAAAAAAAAQM+xEn7bt2/X6NGjz8u4jh8//tzrAAAAAAAAAC4fK+HX1NR0wbeiZmZmKhaLaf/+/T1WMQAAAAAAAAA+K+F3/Phx9e/f/7zjZ79vIuTLLkO+HwFA7zBp0qTLXQUAceJ+Fw2A3uPmm2++3FUAECc8fwNXDivhN3DgQJ04ceK84+3t7ededzHhANE1efLky10FAHFy3XXXXe4qAIgT90czAPQePH8DVw7rZywzMzMv+LHds78y+mmb/8bGRl111VWdjqWlpTHpAAAAAAAAAP9fS0uLWlpaOh07ffp0t8+3En55eXl688031dra2umHO9577z0lJSUpLy/voucPHz68y58ABwAAAAAAAHDhN8i1tbVp165d3Trf+kjvtGnTdOrUKS1duvTcsZMnT2rFihWaMGGChg8f7vw5AAAAAAAAAD3Meoff+PHjVVJSou985zs6cOCAsrOztWLFCu3evVtVVVXxqiMAAAAAAACAbrISfpL085//XHPnztXKlSvV0tKi3Nxcvfbaa5o4cWI86gcAAAAAAADAYCf8+vXrp2effVbPPvtsPOoDAAAAAAAA4BJY3+EHAAAAAAAA4LMtKRaLxeJdyLZt2zRu3DhNnTpVQ4cO7dY5f/krwN1VV1dnxWdmZtpl5OfnW/HHjh2z4kPqVFNTY8VnZGTYZbjtUVtba5fhKioqinsZhYWFVvyiRYvsMkLaw5GTk2Of446l1tZWuwyX296rVq2yy3DvVX19vV1GdXW1Fb9w4UIrfvv27Va85I9vdx6UpKamJiu+ubk57mW44ztkfnbHRki/devljqWQfu5etzsuJP863PiQfu6ufe76LUnZ2dlWvHsdy5cvt+Ilvw+6a4zkj9eQfaR7Tsj66nLXgIqKCruMePeRkPncbW93XEhSQ0ODFe/2D7ftQoTMnW77hfRzt0+VlJTE9e+HnBOyP3fvbWpqqhXvjgvJX8MTsad3x1LIHOIKWTPce+vuQ0Kep91+HrIeu+2RiOeGkDJc7lwYcm+72w8PHjyol19+WVu3blVBQcFFY3mHHwAAAAAAABAhJPwAAAAAAACACCHhBwAAAAAAAEQICT8AAAAAAAAgQkj4AQAAAAAAABFCwg8AAAAAAACIEBJ+AAAAAAAAQISQ8AMAAAAAAAAihIQfAAAAAAAAECEk/AAAAAAAAIAIIeEHAAAAAAAARAgJPwAAAAAAACBCSPgBAAAAAAAAEULCDwAAAAAAAIiQpFgsFot3Idu2bdO4ceP0N3/zNxowYEC3zikpKbHLqa2tteJbW1vtMlJSUqz4wsJCuwxXc3OzFZ+I63brFHKf6uvrrXi3f0j+defk5NhluLKzs6346upquwz3ut14SWpoaLDi3T6SkZFhxUtSXV2dFR/S3u69yszMtOJD+rk734bMIW69mpqa7DLce+W2d8h1V1RUWPGrVq2Kexk1NTV2Ga6Q8edy+5TbP0KuwV2X8vPz7TLc9TVkLLlSU1Ot+KKiIrsMt71D2s9dl9x7W1paasVL0rp16+xzXG4/dOdOd+yFCOlT7nW4fTBkb+uOb3fOkfy9SyLW40SMb/e6Q/bPM2bMsM+JN3fND3kGd+cQt07us48Uth92uf08ZB/pcp8zQvrs9u3brXh3rpX8+TMRz3Bun8rLy7PL6O7YaGtr065du7R161YVFBRcNNZ6h9+WLVtUXl6usWPHKiUlRSNHjtT06dODFhsAAAAAAAAAPa+PE/zss8/q3XffVUlJiXJzc9Xc3KxFixapoKBAmzZt0s033xyvegIAAAAAAADoBivh961vfUv/9V//pT59/nxaaWmpbrnlFn3/+9/Xz372sx6vIAAAAAAAAIDusxJ+EyZMOO9Ydna2xowZo507d/ZYpQAAAAAAAACE6ZFf6T1w4ICGDh3aE38KAAAAAAAAwCW45ITfypUr1djYqHvuuacn6gMAAAAAAADgElxSwu/DDz9UeXm5Jk6cqAceeKCn6gQAAAAAAAAgUHDC78CBA/rKV76itLQ0VVdXKykpqSfrBQAAAAAAACCA9aMdZx09elR33HGHjh49qrffflsZGRndOu/QoUNKTu6cY0xNTVVqampINQAAAAAAAIDIaWlpUUtLS6djp0+f7vb5dsLvxIkTuvPOO9XQ0KDXX39dN9xwQ7fPHTp0qAYMGOAWCQAAAAAAAFwx0tLSlJaW1ulYW1ubdu3a1a3zrYTfmTNnVFpaqk2bNulXv/qVxo8f75wOAAAAAAAAIM6shN83v/lNrV27VsXFxTp06JBefPHFTq/ff//9PVo5AAAAAAAAAB4r4bdjxw4lJSVp7dq1Wrt27Xmvk/ADAAAAAAAALi8r4ffGG2/Eqx4AAAAAAAAAekDyp4cAAAAAAAAA6C2SYrFYLN6FbNu2TePGjdP06dM1bNiwbp3T2tpql1NXV2fF5+fn22U0Nzdb8SUlJVb88uXLrXjJv46Q6164cGHcy3ClpKTEvQxXbW2tfU5OTo4VX1hYaMWH3Kf6+nor3r0GSaqpqbHi3TotWLDAipek6upqKz6kn4eMcUdFRYV9jlunkPZ2+23IvZ0xY4YVX1paasWHjO+MjAwr3l1jQrhzQlNTk12Gex1FRUV2GS63/UL6ubt3ycvLs8vYvn27FZ+dnW3Fh+y/3Pk5EXuEEO7YcPtUyPzs7m0TsX92ywgZS66QeSozM9OKd/tHSD93+1QirjvkOtw+5cZXVlZa8ZK/10lNTbXLcLnzs3ufJL+93TpJfh9paGiw4t01JoS7X5P8Z5mysjK7DFfIXtXlrgEh43XdunVWvNvPQ+6TO9/Gc+07dOiQ1q5dq61bt6qgoOCisbzDDwAAAAAAAIgQEn4AAAAAAABAhJDwAwAAAAAAACKEhB8AAAAAAAAQIST8AAAAAAAAgAgh4QcAAAAAAABECAk/AAAAAAAAIEJI+AEAAAAAAAARQsIPAAAAAAAAiBASfgAAAAAAAECEkPADAAAAAAAAIoSEHwAAAAAAABAhJPwAAAAAAACACCHhBwAAAAAAAERIn8tdga60trba5xQVFcWhJpemurraig+57pSUFCu+rq7OLqOystKKb25utuLr6+uteElqaGiwz3FlZ2db8Tk5OXYZbvvNnj3biq+oqLDiJSkzM9OKD2m/wsJC+xxHSP9w28KNl6Sqqiorvqamxi7DFe+2kPx5x51DpLC+7gi5T24fyc/Pt8uYMWOGFV9SUhLXeMnv5yFzpzs23Parra214kOEzCEZGRlWvDv23LaTEnOvErGviPd4DZnX3LHh7jslacqUKVa8297uXkry67Rq1Sq7DJc7h4SM7wULFljx7vwv+c8aTU1NdhmusrIyKz5kfLvcuVbyx6s7NkL2Ie54DXkWdbl1SkRbhJSRmppqn+MIeWZPRF7AvbchewR3fU3EuuTOhSHPJd1tj927d2vt2rXdir3kd/jNnz9fycnJys3NvdQ/BQAAAAAAAOASXVLCr7GxUd/73veC/gcLAAAAAAAAQM+7pI/0futb39Lf//3f69SpUzp8+HBP1QkAAAAAAABAoOB3+L311ltavXq1XnjhhZ6sDwAAAAAAAIBLEJTwO3PmjCoqKjRr1iyNGTOmp+sEAAAAAAAAIFDQR3oXL16sPXv26Pe//31P1wcAAAAAAADAJbDf4XfkyBE99dRTevLJJ5Wenh6POgEAAAAAAAAIZCf8nnjiCQ0ZMkTl5eXxqA8AAAAAAACAS2B9pLehoUHLli3Tf/7nf6qxsVGSFIvF1N7ero6ODu3evVuDBw9WWlraBc+vra1V//79Ox0bPXq0Ro8eHVh9AAAAAAAAIFree+89bd68udOxtra2bp9vJfwaGxsVi8VUUVGh2bNnn/d6VlaW5syZox/+8IcXPL+wsFDDhg1zigQAAAAAAACuKBMmTNCECRM6Hdu9e7f+/d//vVvnWwm/sWPHas2aNecdf+KJJ9Ta2qqFCxcqKyvL+ZMAAAAAAAAAepCV8BsyZIiKi4vPO/78888rKSlJU6dO7bGKAQAAAAAAAPDZP9rRlaSkpJ76UwAAAAAAAAACWe/w68obb7zRE38GAAAAAAAAwCXqsXf4AQAAAAAAALj8euQdft01aNAgpaSkdCu2qKjI/vu1tbX2Oa7q6morvqSkJK7xklRTU2PFt7a22mU0NDRY8Xl5eVZ8fn6+FR8ipE9VVlZa8QsWLLDLcM/53ve+Z8U3Nzdb8ZLfRzIyMuwy3HPcPrJq1SorXpIqKiqs+KqqKrsMdw7p7px51owZM6x4yZ9DQrh9KmROcPu6e29D1hi3ny9cuNAuIycnx4pPTU214kOu251v6+rq7DLc9nP7hzsfSNLy5cut+MLCQrsMdyy5152ItTJkfLv31u0fkj+W3HtbX19vxUvx339Jfr3cfhvSFiFzoSvee7yQvVF2drYVH7Knd/t5yB7B3Yu4a0BmZqYVL/n9NmT/3NTUZMW74zWkLdzrCHkWdfcJ7n0KWStdIfsQd352ryNkfLvzbchYcoWsfW57uGWErEtue7j7Fqn7feqTTz7p9t/kHX4AAAAAAABAhJDwAwAAAAAAACKEhB8AAAAAAAAQIST8AAAAAAAAgAgh4QcAAAAAAABECAk/AAAAAAAAIEJI+AEAAAAAAAARQsIPAAAAAAAAiBASfgAAAAAAAECEkPADAAAAAAAAIoSEHwAAAAAAABAhJPwAAAAAAACACCHhBwAAAAAAAEQICT8AAAAAAAAgQkj4AQAAAAAAABGSFIvFYu5J27Zt09NPP6133nlH7e3tysrK0j//8z+rvLy8y/hx48Zp6tSpGjp0aLfKyM/Pd6ulpqYm+5x4q62tteILCwvtMtx71draapfhXkdzc7MVn5KSYsVLfnunpqbaZbj1qq+vt8vIzMy04vPy8qz4hoYGK15KzL11ZWdnW/Ehc4jb3u64kPzx584JIX3QrZPbZyWprq7Ois/IyLDLcK/D7SMh99YdfyFrgHuvqqurrXh37En+WKqpqbHLyMnJseKLioqs+ESMJfcaJH/ecdsvZF5z54SQNd9tv8rKSrsMd05wryNkXnPvrTvXSv74q6iosMtwufdq+fLldhmlpaVWvHtvQ/bbrpD52V2X3LEn+X3K3XeGrJXudYfsI1etWmXFl5WVWfEha2VJSYkV7+4RJL89ErHfdueQkL3tunXrrHh3zQ9ZMxIxvt32+Czu6UPGkjtPhayV3V1nmpubVVVVpa1bt6qgoOCisX3cSvz2t79VcXGxCgoK9OSTTyolJUV//OMftW/fPvdPAQAAAAAAAOhhVsLv2LFjevDBBzV16tSg7D8AAAAAAACA+LK+w+/FF1/UwYMHNX/+fElSW1ubAj4RDAAAAAAAACBOrITf66+/rsGDB2vv3r268cYblZKSosGDB+uRRx7RiRMn4lVHAAAAAAAAAN1kJfzq6+vV0dGhu+66S0VFRVq9erVmzpypJUuW2F86CgAAAAAAAKDnWd/h19raquPHj+vhhx/W888/L0m6++67deLECS1dulTPPPOM/vZv/zYuFQUAAAAAAADw6ax3+A0cOFCSdM8993Q6ft999ykWi2njxo09VzMAAAAAAAAANusdftddd50++OADfe5zn+t0fNiwYZKklpaWi56/efNm9evXr9OxrKwsZWVlOdUAAAAAAAAAIusPf/iDPvjgg07HnN/PsBJ+48aN04YNG9TY2KicnJxzx/fv3y9Juvbaay96/vjx4zV06FCnSAAAAAAAAOCKMmbMGI0ZM6bTsebmZlVVVXXrfOsjvaWlpYrFYvrJT37S6fiyZcvUt29fTZ482flzAAAAAAAAAHqY9Q6/vLw8lZWVqaqqSh0dHZo0aZLeeOMN/fKXv1RlZaUyMjLiVU8AAAAAAAAA3WAl/CTpxz/+sUaOHKmqqiq9+uqrGjlypF544QXNnj07HvUDAAAAAAAAYLATfldddZXmzp2ruXPnxqM+AAAAAAAAAC6B9R1+AAAAAAAAAD7bkmKxWCzehWzbtk3jxo3T1KlTu/0rvX/5K8Dd1dTUZMXX1dXZZbhaW1ut+Pz8fLsM95yUlBS7jObmZvscR21trX1OYWGhFe/2D8lvv3jfJ0k6duyYFe/ep5BzFi1aZJfhfufnZ/E7QlNTU+1z3L7u9tuQOcSdE0L6uTuWioqK7DLifW9DxpLbR0LmQvdeudfttp3k95GQ8e3Wy91XhOwR3PEXUkZ2drZ9jiOkD5aVlVnxIfuQ6upqK76iosIuo76+3opPxJq/fPlyKz5knnKvw+2DIffJHUshezy3r7tzbchYmjJlihW/bt06uwy3PRIxltx+W1NTY8VL/rwTMk9lZmZa8e46lpeXZ8VL0oIFC6x4tw+GSMQzgNtrIIzuAAAIyklEQVRHQsar297udSQiHxLCve6FCxfaZVRWVlrxbvuFPC+57RfPnMvevXv1gx/8QFu3blVBQcFFY3mHHwAAAAAAABAhJPwAAAAAAACACCHhBwAAAAAAAEQICT8AAAAAAAAgQkj4AQAAAAAAABFCwg8AAAAAAACIEBJ+AAAAAAAAQISQ8AMAAAAAAAAihIQfAAAAAAAAECEk/AAAAAAAAIAIIeEHAAAAAAAARAgJPwAAAAAAACBCSPgBAAAAAAAAEULCDwAAAAAAAIgQEn4AAAAAAABAhPRxT2hoaNC//du/6Z133tGRI0c0YsQI3XffffrXf/1XDRw48KLnTpw4UVlZWd0qp66uzq2aioqKrPjW1la7jMLCQiu+trbWLsNVX19vxefk5NhluNeRkZER13hJOnbsmBXv9g9Jqq6ujnsZ7rW7bZGfn2/FS1JNTY0VX1JSYpcR77HhjlVJam5utuLd+yRJ2dnZVnxeXp4VH3JfFy1aZMUvX77cLiMlJcWKD7kOt6+7bTFlyhQrXpLWrVtnxc+ePdsuY8GCBVa8e59C1kq3vUPWfLde7loZMnc2NTVZ8SH31pWammrFh+wR3LUyZM13hfQpd053x6t7nySpsrLSiq+oqLDLcOcQlzsfhMjMzLTPcce4OzZC2qKsrMyKLy0ttctw551Vq1bZZbjcPhKy33b3FYmYCz+LYyNk/zVjxgwrfuHChXYZLndshLSFO5YSsRd297YhzzLu3BnyTObeK3eP58ZL/pzuPldK3W/zlpaWbv9NK+G3b98+/d3f/Z3S0tI0e/Zspaena+PGjXrqqae0bds2rVmzxvlzAAAAAAAAAHqYlfD72c9+pqNHj2rjxo268cYbJUlf//rXdfr0af385z/XJ598omuuuSYuFQUAAAAAAADw6azv8Dv7Ecphw4Z1Op6RkaHk5GT169ev52oGAAAAAAAAwGYl/CZPnqxYLKaysjLt2LFD+/bt08svv6wlS5Zozpw5n/odfgAAAAAAAADiy/pI7+2336558+ZpwYIF+tWvfiVJSkpK0hNPPKFnnnkmLhUEAAAAAAAA0H32r/SOGjVKkyZN0rRp05Senq7XXntN8+fPV0ZGhh555JF41BEAAAAAAABAN1kJv5deekkPPfSQGhoazv3E9913363Tp0/r8ccf17333qu0tLS4VBQAAAAAAADAp7MSfosXL1ZBQcG5ZN9ZxcXF+ulPf6q6ujr9wz/8Q5fn//SnP9WgQYM6HZs4caK+9KUvOdUAAAAAAAAAImvPnj3as2dPp2MdHR3dPt9K+B04cEDp6ennHT9b4KlTpy56/oMPPqisrCynSAAAAAAAAOCKMmLECI0YMaLTsZaWFm3YsKFb51u/0jt69GjV1dWpoaGh0/Ff/OIXSk5OVm5urvPnAAAAAAAAAPQw6x1+jz32mH7zm9/oS1/6ksrLyzVkyBCtXbtW69ev16xZs5SRkRGvegIAAAAAAADoBivhV1hYqHfffVdPP/20Fi9erMOHD+v666/XggUL9Nhjj8WrjgAAAAAAAAC6yUr4SdIXvvAF/frXv45HXQAAAAAAAABcIus7/AAAAAAAAAB8ttnv8LsU77zzjj766KNuxaakpNh/v7q62orPycmxy3DrlZ+fb8XX1tZa8ZJ/HSFlHDt2zIp371NRUZEVL8n+zsjKykq7jNLSUis+pN/W1NRY8a2trXYZrubm5rjGS/69+usfC+rpeEnKzs624kPa222/zMxMK76srMyKl6S6ujorPqS9XW5bSFJeXp4V715HyNyZmppqxW/fvt0uw11nCgsLrfiFCxda8ZI0ZcoUK969BslfA9yxtHz5cite8tdjt06SVF9fb8W7/dztH5J/r0L6lLvHa2pqssuYPXu2FV9VVWXFh/Rz91651yD565I7F4bM5+7cGdLeFRUVVrzbFu7aKvnzWkif+s53vhP3MkLOcYTsQ9z2CHlucIXsK+ItZA5ZtGiRFe/eW3f+l/znq5DfInDbz92ju8/fkj823HlQ8q97xowZdhlum7t90J0Hpfjvv6TuP1u2t7d3+2/yDj8AAAAAAAAgQkj4AQAAAAAAABFCwg8AAAAAAACIEBJ+AAAAAAAAQISQ8AMAAAAAAAAihIQfAAAAAAAAECEk/AAAAAAAAIAIIeEHAAAAAAAARAgJPwAAAAAAACBCLnvC73/+538udxUAxMmOHTsudxUAxMn//u//Xu4qAIgT1m8guvbs2XO5qwAgQUj4AYgbHhiA6CLhB0QX6zcQXST8gCvHZU/4AQAAAAAAAOg5JPwAAAAAAACACCHhBwAAAAAAAERIn0QUcvz4cUnSn/70p/NeO3nypA4dOnTe8ba2Nrsc95xBgwbZZSQneznSw4cPW/EtLS1WvCQ1NjbGvYzW1ta4lrF7924rXvLrFNKn3Ht74sQJu4wL9f+LOTueuivkezLd9uvqPrW3t3f5mnvdIf3W1dzcbMW71xDCnadC+qDrs9gWktTQ0GDFu/NzIto7hFsvd07oaq49depUl6+5c2fIeuzO6UePHrXiQ9o75Dpc8e6HiViP//CHP9hluH0qZM1355CeWisvxr2OkLnzQmPjYuv3J598Yv39gwcP2nU6c+aMFR8yLnbu3BnXMhLxLPPRRx/ZZbjjNaT9QuYRR8i9da/7/ffft8twr9vtUyFrzIXmqY6Oji7nL3celPw5wS0jZHy7derTx0+LtLe3W/HuGuD2Wclfl/bu3Rv3MkL6lNvm7lwYcm8/S8+JZ/Nq3ckLJMVisVjcavL/vfjii/rqV78a72IAAAAAAACASFu5cqXuv//+i8YkJOF36NAhrV+/XqNGjdLAgQPjXRwAAAAAAAAQKcePH9fHH3+s22+/XUOHDr1obEISfgAAAAAAAAASgx/tAAAAAAAAACKEhB8AAAAAAAAQIST8AAAAAAAAgAgh4QcAAAAAAABECAk/AAAAAAAAIEJI+AEAAAAAAAARQsIPAAAAAAAAiJD/B4p9O6ICnDJ6AAAAAElFTkSuQmCC", "text/plain": [ "PyPlot.Figure(PyObject )" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "matshow(spikes', aspect=\"auto\", cmap=\"gray\");" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Number of parameters\n", "- $a$: 2$U$\n", "- $\\sigma$: 2\n", "- $\\eta$: 2$TU$" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "-21631.47778954906" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "npars = U * VinDsl.num_pars_advi(RReal()) + VinDsl.num_pars_advi(RPositive()) +\n", "T * U * VinDsl.num_pars_advi(RReal())\n", "xx = 0.1 * randn(npars)\n", "ELBO(xx)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " 1.562912 seconds (8.04 M allocations: 679.587 MB, 6.21% gc time)\n" ] }, { "data": { "text/plain": [ "2022-element Array{Float64,1}:\n", " 96.2983 \n", " 38.2569 \n", " 2615.62 \n", " -7682.52 \n", " -417.951 \n", " -374.225 \n", " -142.068 \n", " -58.1995 \n", " 942.229 \n", " -725.058 \n", " 273.02 \n", " -36.452 \n", " 356.29 \n", " ⋮ \n", " 17.0867 \n", " -16.259 \n", " 20.1646 \n", " -3.70256\n", " 4.73903\n", " 3.36335\n", " 5.88467\n", " 7.44689\n", " 13.5717 \n", " -6.19642\n", " 8.93425\n", " -12.4798 " ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "∇L = (storage, x) -> ForwardDiff.gradient!(storage, ELBO, x)\n", "stor = similar(xx)\n", "∇L(stor, xx)\n", "@time ∇L(stor, xx)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "yy = copy(xx)\n", "gg = similar(yy)\n", "avg_sq_grad = similar(gg)\n", "firstpass = true\n", "elbo = Float64[];" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "collapsed": false, "slideshow": { "slide_type": "-" } }, "outputs": [ { "data": { "text/plain": [ "\"Iteration 1 -21901.643192506774\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 2 -21784.408566088663\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 3 -22761.50851929847\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 4 -20361.289769864685\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 5 -21527.27940310702\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 6 -17984.190834400826\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 7 -17219.647018742507\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 8 -20096.171998652517\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 9 -15719.133230555724\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 10 -14807.092746386164\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 11 -14563.375832332598\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 12 -14043.741827632271\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 13 -13715.407754031177\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 14 -12958.378054523191\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 15 -11893.664828330544\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 16 -11416.167999849515\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 17 -10754.339545075985\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 18 -10049.352149340006\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 19 -9402.409352119079\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 20 -9658.607301884787\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 21 -8327.828020622168\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 22 -8908.349523085642\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 23 -7730.0682606714045\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 24 -6830.5324054293405\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 25 -7717.520227847391\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 26 -6949.580812127499\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 27 -6652.553936474353\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 28 -5976.654400816284\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 29 -6357.3658596257765\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 30 -5812.548575407791\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 31 -5490.162040191422\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 32 -5059.303556650625\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 33 -4781.426677334765\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 34 -4551.261856632452\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 35 -4308.270095528986\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 36 -4292.079328215761\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 37 -3994.083445482836\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 38 -4023.492487795364\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 39 -3941.4627487592675\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 40 -3773.625767885166\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 41 -3554.285420824762\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 42 -3471.2177896320422\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 43 -3606.801729149118\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 44 -3458.755837296366\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 45 -3312.558257314446\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 46 -3422.1017210636055\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 47 -3486.693965133956\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 48 -3212.4499594516324\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 49 -3295.4129262072634\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 50 -2894.734068845817\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 51 -3044.5778606475983\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 52 -3043.3890236880393\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 53 -3096.705300079607\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 54 -2800.2168429386597\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 55 -2823.6872893306668\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 56 -2724.815827991035\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 57 -2536.2819942592146\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 58 -2813.8252996675283\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 59 -2622.910766909394\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 60 -2545.5242115244478\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 61 -2716.8188182666395\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 62 -2547.8211304211463\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 63 -2487.1965455038285\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 64 -2595.838009852493\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 65 -2555.6287417482945\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 66 -2496.6940447238594\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 67 -2440.753239815258\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 68 -2529.6677493249813\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 69 -2342.0300345334454\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 70 -2422.372643545843\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 71 -2468.877089588612\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 72 -2275.624680316577\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 73 -2340.410675965846\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 74 -2363.202887175471\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 75 -2298.6586223903637\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 76 -2424.8088286744824\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 77 -2254.9648424563484\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 78 -2256.767659962991\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 79 -2268.489729925305\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 80 -2321.730033369412\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 81 -2263.3048735619204\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 82 -2247.089332315065\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 83 -2239.594407428411\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 84 -2157.729730733962\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 85 -2320.5130380270866\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 86 -2181.973441053312\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 87 -2151.7475185448898\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 88 -2136.0583167507343\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 89 -2207.967863877849\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 90 -2213.9749887090247\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 91 -2239.828309206509\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 92 -2176.8503933645948\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 93 -2127.587617212332\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 94 -2206.123110053107\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 95 -2096.1898543678385\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 96 -2110.604683041632\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 97 -2142.6319010119396\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 98 -2072.0808398411173\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 99 -2082.9017019162893\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 100 -2067.3815735554995\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 101 -2038.6749828045126\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 102 -2008.7349778361252\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 103 -2175.2734896075244\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 104 -2163.3208904355224\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 105 -2062.4252397849423\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 106 -2027.6384133912975\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 107 -2028.4687089222305\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 108 -2092.8963261403974\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 109 -2044.4750625303006\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 110 -1980.050522057777\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 111 -2059.126407304515\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 112 -2016.3900277082657\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 113 -1961.0534237732666\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 114 -2076.216579276364\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 115 -1965.2993420464118\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 116 -2048.9245640823733\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 117 -2052.929140493427\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 118 -2007.0027668345772\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 119 -1930.4442542333984\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 120 -2060.5535192709685\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 121 -1949.518181457218\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 122 -1952.78904923031\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 123 -2000.9738932013327\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 124 -2099.798320407929\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 125 -2093.411012505818\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 126 -1966.694448941153\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 127 -1943.7472773324391\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 128 -2001.3031484555793\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 129 -2032.0800170265766\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 130 -1923.241793347899\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 131 -1910.3122495794478\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 132 -1996.2578334156574\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 133 -1957.6409762797805\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 134 -1951.2053297682116\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 135 -1976.5174852613877\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 136 -1908.6343265168282\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 137 -1925.3980986570036\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 138 -1887.3830020508765\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 139 -1918.5944250428877\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 140 -1991.5230883643699\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 141 -2025.4096647444849\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 142 -1925.7028873309512\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 143 -1950.1036384208073\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 144 -2057.729579572335\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 145 -2017.467800560341\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 146 -1859.708833399075\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 147 -1888.4565974368297\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 148 -1939.0910823624286\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 149 -1948.19136384097\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 150 -1919.7321499975917\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 151 -1950.7692580396708\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 152 -1890.3043823100202\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 153 -1905.7101972997893\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 154 -1888.98786569496\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 155 -1905.1120590604025\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 156 -1940.6845582157573\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 157 -2026.209087976532\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 158 -1931.0174026672742\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 159 -1857.8731177192637\"" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "\"Iteration 160 -1905.7424392670107\"" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "step_size = 0.1\n", "decay = 0.9\n", "eps = 1e-8\n", "nsteps = 160\n", "\n", "for jj in 1:nsteps\n", " ∇L(gg, yy)\n", " nn = norm(gg)\n", " if nn/length(gg) > 10.\n", " gg /= norm(gg) # crude gradient clipping\n", " end\n", " if firstpass == false\n", " avg_sq_grad = avg_sq_grad * decay + gg.^2 * (1 - decay)\n", " else\n", " avg_sq_grad = gg.^2\n", " firstpass = false\n", " end\n", " yy += step_size * gg ./(sqrt(avg_sq_grad) + eps)\n", " push!(elbo, ELBO(yy))\n", " display(\"Iteration $(length(elbo)) $(elbo[end])\")\n", "end" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false, "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAtgAAAIUCAYAAAAg+VvKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3X1Y1FX+//HXwCgCCt6VIalpea+F4N2qGd5l2m5kqamVe+2Wml7bjX1rK7NyW/OXabrbpttmlpVpW1ZqmbVRma7iDSh5E1KalSi2uRoCKQh8fn+cnYGRQUGBz2fg+biuuWY5c+YzZz55Xfvi8D7nuCzLsgQAAACgUgTZPQAAAACgJiFgAwAAAJWIgA0AAABUIgI2AAAAUIkI2AAAAEAlImADAAAAlYiADQAAAFQiAjYAAABQiQjYAAAAQCUiYAMAAACViIDtR35+vh566CFFR0crLCxMvXv3VmJiot3DAgAAQAAgYPvx29/+Vn/5y190++2367nnnpPb7dbw4cO1adMmu4cGAAAAh3NZlmXZPQgn2bp1q3r37q1nn31WU6dOlSTl5eWpS5cuatasmf7973/bPEIAAAA4GTPYZ1ixYoXcbrcmTJjgbQsJCdEdd9yhpKQkHTp0yMbRAQAAwOkI2GdITU1Vu3btVL9+fZ/2nj17el8HAAAAykLAPkNmZqaioqJKtUdFRcmyLB0+fNiGUQEAACBQuO0egNOcPHlSISEhpdrr1avnfd2fo0eP6uOPP9Zll12m0NDQKh0jAAAAKu7kyZP67rvvNHToUDVt2rTKPoeAfYbQ0FDl5eWVaj916pT3dX8+/vhj3XbbbVU6NgAAAFy4pUuX6tZbb62y6xOwzxAVFeW3DCQzM1OS1Lx5c7/vu+yyyySZ/2AdO3assvHVNFOnTtX8+fPtHkbA4b5VHPfs/HDfKo57dn64bxXHPau4tLQ03Xbbbd7cVlUI2GeIiYnRunXrlJOT47PQcfPmzXK5XIqJifH7Ps/MdseOHRUbG1stY60JIiMjuV/ngftWcdyz88N9qzju2fnhvlUc9+z8VXU5L4sczzBy5EgVFBToxRdf9Lbl5+dryZIl6t27t6Kjo20cHQAAAJyOGewz9OzZU6NGjdIjjzyiH3/8UVdccYWWLFmi77//Xq+88ordwwMAAIDDEbD9eP311/XYY49p6dKlOn78uK688kqtWbNGffv2tXtoAAAAcDgCth9169bV7NmzNXv2bLuHUuONHTvW7iEEJO5bxXHPzg/3reK4Z+eH+1Zx3DPnclmWZdk9iJpg+/btiouLU0pKCgsOAAAAHKi68hqLHAEAAIBKRMAGAAAAKhEBGwAAAKhEBGwAAACgEhGwAQAAgEpEwAYAAAAqEQEbAAAAqEQEbAAAAKASEbABAACASkTABgAAACoRARsAAACoRARsAAAAoBIRsAEAAIBKRMAGAAAAKhEBGwAAAKhEBGwAAACgEhGwAQAAgEpEwAYAAAAqEQEbAAAAqEQEbAAAAKASEbABAACASkTABgAAACoRARsAAACoRARsAAAAoBIRsAEAAIBKRMAGAAAAKhEBGwAAAKhEbrsHAAAAAMmyJJer+j/38GHp5Zeliy+WrrhCattWio6WggJsGnbvXuntt6WICKlfP+mqqyS3TUmXgA0AAKrU0aMm9NSta/dIzu7kSemjj6SffpJ+/tk86taVbr9duvzy879uYaH0739LH3xgfr7oIhNmGzSQ0tOl1FTzOHBAatVK6tBBat/eBMQbbpAaNjz3Z3z7rRl3q1bm2uUNx+vXS6NHS9nZUl6eGask1a8v9e0rDRggxcdLsbFSnTq+783Kkvbtk/bvN6916GDuU9265rWPP5bef19KTDS/ODRpIjVtKkVGmtd/+sk8srKkRo3MfWnaVIqKkrp0kbp2NY+WLcv+xePoUWnZMun116XkZHNP8/PNd6lfX+rTR3r4YfM9qhMBGwCAGuj776XwcBNYzua//5U+/VTavl0aN0668srSfb78UtqwQbrpJql58/J9/i+/SO+8I73yivT55yYwvvSSmVn0sCzprbekv/7VhFu32wS1sDCpZ0/p6qtNQIqMLN9nWpa0cqX0xRfSwYPmceiQ1KaNNGyYecTElA5rhYXSq69KTzwhZWSY1yMjTeg7flx68knp+uulu++WrrlG+uoraccO8zh0SDp1qvjRoIF02WXmceml0tat0ooV0pEjZlY4NLQ4VErmc2JipOHDzTi//96E7lWrpPnzTVi94Qbpttuk664r/UvKnj3SzJnSP/9pvr8khYSYUHrTTdJDD5nv4e9e/eUv0oMPmvv8z3+aIP/ddyY0f/mluY9//rMJqJJUr54JreHhUm6uCbdnCg42n33woFRQYP49/fa35nsfPWr+vf38s3TJJSY8X3SR+eXLE7iPHjX3YM2a4nvUvLk0aZI0caJ5n2TGOXeutHix+e83fLgZ5/XXm/9+ycnm3+yGDcX3pTq5LMuOj615tm/frri4OKWkpCg2Ntbu4QAAbJCdLX3zTfFj3z6pf3/pd7+r2J/+8/Olr782M4OFhea9LpcJQH36+J8Jtixp1y7pvfekd9+Vdu40YTUhQbrzTmnwYBN+jhwxoW/zZjOzmJxs3hsZacY/YYIJlBdfbELSY49Jr71m+rjd0ogR0pQpJmiW/E6FhdLu3dKmTWa29v33zfXi46VbbjHXSEoy7/1//8/M2D7wgLRtm3TttWbms6DAPI4fN31//NHMxPboYfredFPZM7Pffmuu/fHHUrt2Zia3ZUszG7p7t/muOTkmoHXrZgJ/+/YmzM+ebULz6NHSn/5k3u/5nJMnzQzp3/5mQqfLVVzK0a6dCcWhoSZ8hoSYUPjdd+Zx7JgJ2aNGmXvQs2fxPcvLk06cML8AlfVv4/BhaflyMzv75ZfmMzp2lDp3No/t2014b9FCmjZN6tVL+uEHE1D37jW/NNSpY4Ln3Xebcf7wg5SSIi1dav6tPPigNGtW2aUUp0+b/l99ZUJ1To551KtXXE5y+eXm3+zeveaXg6+/Nvfl1782/w3Oh2WZX3Z27ZJWrzb34PRpcy+Dgsx9adhQuvdeafLkc/8i6VFdeY2AXUkI2ABQs1iWCamff26CU0SEmZ2sW9cEp6NHzePQoeJAfeRI8fubNDHh6ssvzezb3/5W+k/sJW3dKj3/vAlN6ekmaPoTGWmCy003mYCzaZP5M/+GDSaQRESY12+8UcrMlBYtMgGzRQsT5H74wVzn4ovNn82HDpWGDDE///3v0owZUlGRCeZvv22+84wZJnwuXy4tXGiCVOPGJpzWqWMemZkmULvdJsAOG2ZmLtu0MZ9XWCgtWGCCoNttgmhcnJmFjI/3f//37TPfa/lyE5A7d5YefdSMJTjYXPOXX8yY/vQnMxu6YIH5/mfKzzfB/5NPzKxvenrxLzADBpiQ3aPH2f89bNxo3tu1q5mZrV+/7P6SCaJhYZVTy7xrl7kHX31lxrBnjwmVjzwijR/v/5euH380M9D/+IeZxS4qMjPIkvm3+Ze/SDfffOFjqw7Hj0tLlpj/vqdPS/ffb35xDA+v2HUI2AGGgA0AFVNUZP5P0/Mn486dK/5/lmc6ccKEoDp1TGgKDi5+raDA/B/0n/5kAnLjxubRpIkJgR07mscll0gffmhC3VdfmXBdWFg68NarZwLOJZeYWbx27cyz5+H5s/zLL0t33WVmnlesKD3TtnWrGdOHH5prDB5s6k+7dDHXqVPHhDvPjN6qVWaGetcu8/7gYBNU+/eXBg0y3zskpPj6lmVmiV9/3Yy5Z08z0+kJ3Gf6739NoH7nHemOO8wMZ0SE7/XWrTPBPj/fhJ3Tp8197NNH6t7dhMqyfP+99OyzUu/e0pgx5Q+fSUkmLK5da75HYaH5XMlc4777zH08V+gtKT/flCU0b27P4sIL4Ulv5Rn3t9+aX/AaNjT/VuLizMx+bUTADjAEbAC1VUqKmSE7eVKaN8/MIp4pNdXMvn39dfGfkP/zHxOyPcLDTfnBbbeZoFie1f/5+WaG+eOPTb1oamrxNVu0kH7/e1OesXOnqUVNS5PGjjUh89gx8/jpJzNTunevmQ0tOZaxY83srttt6ms9C8GaNDl7iDzTv/9tZpzr1zezr1lZ5nHwoHmtQwfp8ceLZ2bLY98+E7i7d69YqAx0ycnmnoWEmKBdr55ZDNipk90jQyAgYAcYAjaAmuzkSTOz+csv5nHypCk7ePFFE3guvdQEz9BQs8isWzfzvrw8M6s4e7Z5rV07U/fatq15T9Om5hEWZnZvWLrUBPAmTaRmzYrrWsPDTS1n69bmIZkdGT780Mxat2hhygyuucbM5B4/bhY/LVtm/kwvSQMHSs88Y2bv/CkqMoH1hx/MjgkVCdDl8f335k/5GRmmzMOziG7UqIoFawDnr7ryGruIAIBDnD5tAuvFF5v61XPNShYVFYfaM2Vnm9nbjz4yfz73PEJDTbCLiDB/Lu7a1ZQk9OtnAmVBgbRliwmun31m6mp/+ql4Zrckl8uMc/Vq85yZaep++/Y14bZzZ7O9WVqa2eXgwQfPPivdvbupr92+3SyQy8oys8Z5eeb7eBY7eXYuiI2V/u//zGd27Vr6T+U9e5pShFWrTFgfNOjsf04PCjIh/nwXZZ1Lq1Zmlh1AzUfABoAq9OOPJii2a3f2fnl5ph511SpTW1mvnll8lpBg/vTdooUJiZZlFrS9+66ZKf7pJ7No7OGHi+tuf/hB+s1vzJ66Dz9sQnVRUfGCsBMnzJiOHTMzxnPmmAVSsbFm9vj4cTOrPGSImfX17E3bpImZSQ4LM9f07OXr0aKF+dP9xIlmuze325Q+bN1qtiErD5eruEa0LDk55nuU/Oyy1K8v3Xpr+T4bACoLARsAqsiePSakZmaaGeKJE6WRI0vPOP/yi1nJ//nnZoa2U6firdZ+//vifm63CdG5uaa84uabTTCeOdMsyHvhBXPthATzvGmTWSh3NpZlZpg//dQsDrz2WrOfbPfu51eyEBpqtmP71a9M+C8Z/CtL/fq1q+YYQOChBruSUIMN1D6WZbZgO3jQBNOSQXLLFhNUW7QwZQyvvWYW+TVsaGaX+/Y1j5YtTSDeutWE60GDfD/j+HFTu5uRYT4nK8uUdMTFFZc77NljtoHz7J7Ro4cJ6OWZ4QWA2oQabACoBvv3m31xu3UzW5edy/Hjpv/ataa+2bPv8SWXmEMu7rrLhO4bbzRlER98YEL17bebz1q8WPrXv8ziu8JCMysdFmba+vYt/XmNGpnH2UosOnc2ZSMvv2x253jySVNiAgCwBwEbgOPl5JjFf7//vQmr5/LVV2ZxXFmTE3v2mBnl99835REeI0aYE806dChuO33aBOaPPjKPpCRTz9y5s9lObtgwU4u8cKE5nW7WLPP6wIFmH+GSO1Fcfrl5fdYs8522bjVb3F13nVmkdyGCgsyhCwAA+xGwATjasWOm1GLLFlNT/Pbb/neCyMqS/vlPM4u7ZYtpmzXL1ACX7L9smQnqkZHS9ddLTz1lyjJWr5amTzfBefx4MwOckmLCdX6+Oc1u8GBT53zddab0o6S//93UQr/4otnObtYs/yeredSvb0L4wIEXfo8AAM5CwAbgWJmZprY5M9OUPTz+uDn2eeJE334vvSTdc4/ZiWPYMDNzvHOn2V0jPd2EXrdbeuIJE4LHjzdtJWumb7vN7Ef897+brd0iI02d8223mQV/3bufPTBLZpeNRx6p/PsAAAgsBGwAjnTggJkxzsuTNmwwR1gfPizde6+pVe7c2fR79lnpgQdMecSf/mSOPJbMqXlt25rZ6gMHTBnHu++aA08efND/LHhIiDlu+b77qu97AgBqniC7BwAAJf3yi/SXv5gFh0FBZl/ljh3Na/PmmTrmMWPMSYJPPGHC9bRpZkbaE649br3VHJby1VdmUeK770p//OPZDxsBAOBCMYMNwBabNpka55YtTT3zxRebvZznzjU1zLffLj39tDlcxSM0VHrzTbMNXbdupvzj6aelhx4q+3P69jUnAJ48WXzENgAAVYmADaDapaebhYX5+WbHDY86daTf/c4sTCwrDHfpIs2fb7bEW7DAPJ/LJZdUzrgBACgPAjaA83LkiNmW7v77pVatyv++ggLpt781s9Y7dphTCX/4wRykEhdXencOf+66y5SJlGfLPgAAqhs12ADOy4wZ0nPPmVKNVavK/745c6Rt26RXX5XCw01pSPfu5mCW8oRrD8I1AMCpCNgAKuz7781+0w8/LF1zjQnH991ndvw4m507zcLEP/5R+tWvqmesAABUN0pEAFTYrFlmn+hHHzWz0M8/b3bzeO89sxd0UZFkWVLjxtKAAabeOibG7D/dvr2Z/QYAoKYiYAOokO++M7PXs2aZ0wgl6e67zW4dL79swrXLZR6HDpmt9Z54wixgtCxzPHjJA14AAKhpCNgAKuSpp8zM9Jm7d8TGmseZCgvNdnyffmr2sO7WrXrGCQCAXQjYAPwqKjKzzx07StddJwUHmxMRlywxe0+Hh5fvOsHBUs+e5gEAQG1AwAbg19y5xQe4tGolTZpkFik2aSJNnmzv2AAAcDICNlBLWJa0eLHZv/ree6UGDcruu3mzWcD40EPSTTdJf/+79OST0qlTZlY7LKz6xg0AQKBhmz6gFjh2zATlCRNMUL78cnMK4unTpfv+/LM5xKVHD+nPfzalHa+8YhYsvv229Ic/VP/4AQAIJARsoIZbv1666irpiy/MNnr790vXX292/ujUSXrhBenHH01fy5LuvFPKypKWLTM7f3g0biyNHOnbBgAASiNgAwGqqEj66COzBV5mZunXCwrMftMDBkitW0tffll8WuIrr5ifO3UyM9JRUdLVV0t33CG9844pJbnssur+RgAA1AwEbCDAHD1qjhtv104aNszs6NG5s/TGG2YGWpIOHpQGDjQlHo8/Ln32WeljyLt2NUecHzliAnXDhuYad99tykkAAMD5YZEjEEA+/1z69a/N7PTo0dLrr5ugfc890m23SStWSCNGmGPLw8OldevMzPTZNG0q/e535pGXJ9WtWy1fBQCAGouADQSIlBQpIcGcmLhsmQnGHm+8Id18s3TXXdLKlaYUZPFiUzddEZywCADAhSNgAwHg669NOUinTtK77xYfUV7STTdJ/ftLycnS0KHmqHIAAFD9CNiAw2VkSEOGmBnrNWv8h2uPpk3NqYsAAMA+LHIEHOz776VrrzX/+1//MqcoAgAAZyNgAw6VlGQOeTl50oTrSy+1e0QAAKA8CNiAA73+uhQfL7VvL23dap4BAEBgIGADDnLwoNlyb/x46dZbpcRE6aKL7B4VAACoCBY5AtXs88/N/tRXXCG1bWuet241R5avWSOFhUnz5pm9rNkJBACAwEPABqrRoUNmO72CAiknx/e1mBjp73+Xxo6VGjSwZ3wAAODCEbCBamJZ0p13SqGh0u7d5sTEffukb76RWrWSevRgxhoAgJqAgA1Uk8WLpY8+MmUgnhMWY2LMAwAA1BwscgSqwXffSVOnSnfcIQ0fbvdoAABAVSJgA1WsqEj6/e/NrPW8eXaPBgAAVDVKRIBKdvCg9PzzUmamdPSoWdi4c6fZci8iwu7RAQCAqkbABirRiRPSdddJR45InTpJTZtKvXpJf/yjNGiQ3aMDAADVgYANVJLCQmncOCkjQ9q8WerY0e4RAQAAOxCwgUry8MPS2rVmlxDCNQAAtRcBG6gES5ZIc+dK8+ebEhEAAFB7sYsIcIE2bJAmTTKHyNx7r92jAQAAdiNgAxdg714pIUHq21dasICTGAEAgMMC9quvvqqgoKBSj+DgYP3nP/8p1X/Tpk3q16+fwsPDFRUVpXvvvVe5ubml+lmWpWeeeUZt2rRRaGiorrrqKr355pt+x3D48GGNHj1ajRo1UmRkpG688UYdOHCg0r8rAt+RI9KwYVLz5tK775qjzwEAABxXg+1yufTnP/9Zl112mU97w4YNfX5OTU3V4MGD1alTJ82fP18ZGRmaM2eO9u3bpzVr1vj0nTZtmmbPnq1Jkyape/fuWrVqlcaNG6egoCCNHj3a2y83N1fx8fHKzs7W9OnT5Xa7NW/ePMXHxys1NVWNGjWqsu+NwJKbK/3611J+vvTFF9IZ/zwBAEAt5riALUnXXXedYmNjz9pn2rRpaty4sb744guFh4dLklq1aqWJEycqMTFRgwcPlmRmpOfNm6e7775bf/3rXyVJd9xxh6655ho9+OCDGjVqlFz/+7v+ggULtH//fm3bts37+dddd526dOmiZ599VjNnzqyqr4wAUlgojRkjpaeb+uuWLe0eEQAAcBJHlYiUlJOTo6KiIr+vZWdnKzExUbfffrs3XEvS+PHjFR4errfeesvbtnLlShUUFGjy5Mk+15g8ebIyMjKUlJTkbXvnnXfUo0cPn3Dfvn17DRo0yOeaqN1mzZI+/FB6+20pJsbu0QAAAKdxXMC2LEvx8fGKiIhQWFiYEhIStG/fPp8+u3btUkFBgeLi4nza69Spo5iYGO3YscPblpqaqvDwcHXo0MGnb8+ePWVZlrevZVnauXOnunfvXmpMPXv21P79+/3Wd6N22bhRmjFDmj6d7fgAAIB/jioRCQsL0+9+9zsNGDBAERERSklJ0bPPPqu+fftq+/btio6OliRlZmbK5XIpKiqq1DWioqL073//2/tzZmammjVr5refZEpIJOnYsWPKy8sr85qevm3btr3wL4qAdPy4OanxV7+SHnvM7tEAAACnqrKAbVmW8vPzy9U3JCREkjRq1CiNGjXK237DDTfo2muvVf/+/fXUU09p4cKFkqSTJ0/6vK+kevXqeV/39C2rX8lrneuaJfug9rEsacIE6cQJ6Y03JLejfjUFAABOUmUxYf369RowYMA5+7lcLqWlpaldu3Z+X+/bt6969eqlxMREb1toaKgkKS8vr1T/U6dOeV/39C2rX8lrneuaJfuczdSpUxUZGenTNnbsWI0dO/ac74VzvfSS9M47pu66VSu7RwMAAM5l+fLlWr58uU9bVlZWtXx2lQXsDh06aMmSJeXq668so6QWLVro66+/9ulvWZYyMzNL9c3MzFTz5s19+q5bt85vP0nevo0bN1ZISEiZ1yzZ92zmz59/zh1QEFh+/FG67z4zgz1ypN2jAQAA5eFvgnP79u2l1vBVhSoL2M2aNdP48eMr5VrffvutLrroIu/PXbp0kdvtVnJyskaWSDynT59WamqqbrnlFm9bTEyMFi9erL179/osdNy8ebNcLpdi/rcNhMvlUteuXZWcnFzq87ds2aI2bdr47FiC2uO556SgIGn2bLtHAgAAAoGjdhE5evRoqbYPP/xQKSkpGjZsmLctIiJCgwcP1tKlS3129njttdeUm5vrc3hMQkKC3G63t37b44UXXlB0dLT69OnjbRs5cqS2bdum7du3e9vS09P12Wef+VwTtceJE+YI9EmTJM4ZAgAA5eGopVp9+vRRt27d1L17d0VGRiolJUWvvPKKWrVqpUceecSn71NPPaW+ffuqf//+mjhxog4ePKh58+Zp6NChGjJkiLdfdHS07rvvPs2dO1f5+fnq0aOH3nvvPW3cuFHLli3zHjIjSVOmTNGiRYs0fPhwPfDAA3K73Zo/f76ioqJ0//33V9t9gHO8+KL0yy/S1Kl2jwQAAAQKRwXsMWPGaM2aNfrkk0/0yy+/KCoqSpMmTdLjjz/uUyIiSd26dVNiYqIeeugh3X///WrQoIEmTJigWbNmlbru7Nmz1bhxY/3jH//Qq6++qrZt2+qNN97wKSWRpPr16+uLL77Q1KlT9dRTT6moqEgDBgzQvHnz1KRJkyr97nCevDxp3jxp/HjpfztEAgAAnJPLsizL7kHUBJ6i+ZSUFBY51hAvvSRNnCilpUnt29s9GgAAcKGqK685qgYbcIrCQumZZ6SbbiJcAwCAinFUiQjgFO+9J33zjTlUBgAAoCKYwQbOUFgozZolDRok9ehh92gAAECgYQYbOMMLL0g7dkgbN9o9EgAAEIiYwQZKOHxYeuQRs+91iS3SAQAAyo2ADZRwzz1SWJj0//6f3SMBAACBihIR4H/ef1965x1p+XJObQQAAOePGWxAUk6O9Ic/SNddJ51x/hAAAECFELABSU8+Kf30k7RwoeRy2T0aAAAQyAjYqPV++cXsHDJ1qtS6td2jAQAAgY6AjVpv5UopO1v6/e/tHgkAAKgJCNioFfbvN4sY/Xn1ValfP+nyy6t3TAAAoGYiYKNWmDFDGjHCHH9e0qFDUmKi9Nvf2jIsAABQAxGwUeMVFEgffmiOQH/sMd/Xli6V6taVRo2yZ2wAAKDmIWCjxtu8WTp2TLrzTumf/zTHoEuSZZnykBtvlCIj7R0jAACoOQjYqPHef1+6+GJpwQKpfXtp2jTTnpIipaVRHgIAACoXARs13gcfSNdfb0pBZs6UPvpI+uILM3sdFSUNGWL3CAEAQE1CwEaN9u230ldfSb/+tfn55puluDjpoYfMkei33SYFB9s7RgAAULO47R4AUJU++MDMXHtmqV0u6emni3+mPAQAAFQ2AjZqtA8+kOLjpQYNitsGDzYB+8QJqXNn24YGAABqKAI2aqwTJ6R166Rnny392nvvme37AAAAKhsBGzXWJ59Ip08X11+XFB5e/eMBAAC1A4scUWN98IEpAWnd2u6RAACA2oSAjRqpsFBas8b/7DUAAEBVImCjRtq6VfrpJ+k3v7F7JAAAoLYhYKPGKSqSHn7YlIb07m33aAAAQG3DIkfUOM89J61fb3YQ4RAZAABQ3ZjBRsDIz5cs6+x99u6VHnlEuvde6ZprqmdcAAAAJRGwERBOnpSiosz+1WUpKDAnM7ZsKc2aVX1jAwAAKImAjYCQnCwdOyZ99FHZfZ55xvR79VUpLKz6xgYAAFASARsBISnJPG/a5P/1b76RZsyQHnqIhY0AAMBeBGwEhE2bpKAgac8e6eefS7++cqXkdkuPP179YwMAACiJgA3Hsywzgz1qlPl58+bSfT79VOrfX6pXr3rHBgAAcCYCNhzv22+l//xHuv12qWlTaeNG39fz86UNG6SrmQIzAAAgAElEQVRBg+wZHwAAQEkEbDiep/76V7+S+vQpXYe9ebP0yy8EbAAA4AwEbDheUpLUoYPUuLEJ2Fu2mC35PBITzWsxMfaNEQAAwIOADcfbtMnMXktS375Sbq60a1fx659+Kg0caBZBAgAA2I1IAkfLyZF27jQz15IUFyfVqVNch52dLW3dSnkIAABwDgI2HG3rVqmoqHgGOzRUio0trsNev96UixCwAQCAUxCw4WhJSVJkpNSxY3Fb377FATsx0RyNfsUV9owPAADgTARsONqmTeZkxpL11X36SN9/Lx06ZOqvBw2SXC77xggAAFASARuOZVlmCz5PeYiHpx575Uqz2JHyEAAA4CQEbDjW119Lx44VB2qPqCipdWtp9mzz88CB1T82AACAshCw4VibNpnSj169Sr/Wp4908KDUubMJ3AAAAE5BwIZjJSWZAB0RUfo1z6w25SEAAMBpCNhwrKSk0vXXHldfbZ6vvbb6xgMAAFAeBGw4Ul6elJYmde/u//WuXaXkZGn48OodFwAAwLm47R4A4E96ulRYaEpEyhIXV33jAQAAKC9msOFIu3eb57MFbAAAACciYMORdu+WLr1UatjQ7pEAAABUDAEbjrR7t9Sli92jAAAAqDgCNhyJgA0AAAIVARuOk5srHThAwAYAAIGJgA3H+eor80zABgAAgYiADcfZvdsckd6xo90jAQAAqDgCNhxn927p8sulsDC7RwIAAFBxBGzY5tAhaf/+0u27d7P/NQAACFwEbNjm/vul66+XLMu3nR1EAABAICNgwzZ79pgj0XfsKG47flw6fJiADQAAAhcBG7YoLJT27TP/e/ny4vY9e8wzARsAAAQqAjZs8f33Ul6e1K6d9OabUlGRad+9W3K7TTsAAEAgImDDFl9/bZ4fe0zKyJA2bjQ/794ttW8v1a1r39gAAAAuBAEbtkhPl0JCpDFjpEsvLS4TYYEjAAAIdARs2CI9XWrb1pSDjBkjvf22dPo0ARsAAAQ+AjZskZ5uSkEkaexY6ehRadky6b//ZQ9sAAAQ2AjYsEXJgN2tm1nUOGOG+ZkZbAAAEMgI2Kh2OTnmFEdPwHa5zCz2d99J9epJbdrYOjwAAIALQsBGtfvmG/PsCdiSCdiS1KmTFBxc/WMCAACoLG67B4DaJz3dPJfc67p9e6l3bykmxp4xAQAAVBYCNqpderp00UVSo0a+7f/6l1Snjj1jAgAAqCwEbFS7kgscS2rQoPrHAgAAUNmowUa1KytgAwAA1AQEbFQryzLHpBOwAQBATUXARrXKzDTb9BGwAQBATUXARrXyt4MIAABATULARrVKTzf7XHOYDAAAqKkI2KhW6ekmXNeta/dIAAAAqgYBG1UmPV267z7p1KniNhY4AgCAmo6AjSrz17+ax913F7exRR8AAKjpCNioEkVF0qpVUufO0ksvSYsWSXl50oEDBGwAAFCzcZIjqkRysnT4sLR8uXn84Q/mGPSiInYQAQAANRsBG1Vi5UqpaVOpTx+pVy8pNVW6807zGjPYAACgJquWEpEjR47o4Ycf1sCBAxUREaGgoCCtX7++zP6bNm1Sv379FB4erqioKN17773Kzc0t1c+yLD3zzDNq06aNQkNDddVVV+nNN9/0e83Dhw9r9OjRatSokSIjI3XjjTfqwIEDfvuuXr1acXFxCg0NVatWrTRjxgwVFhae35evpVaulH7zG8ntlkJCpBUrpCZNpIgIqVkzu0cHAABQdaolYKenp2vOnDk6fPiwrrzySrlcrjL7pqamavDgwTp16pTmz5+vCRMm6MUXX9To0aNL9Z02bZoefvhhDR06VM8//7xatWqlcePG6a233vLpl5ubq/j4eG3YsEHTp0/Xk08+qR07dig+Pl7Hjx/36bt27VqNGDFCjRs31vPPP68RI0Zo5syZuueeeyrnZtQC6elSWpp0443FbdHR0tq10nPPSWf5zw8AABD4rGqQk5NjHT9+3LIsy1qxYoUVFBRkffHFF377Dhs2zIqOjrZycnK8bS+99JIVFBRkffLJJ962Q4cOWXXr1rXuuecen/f379/fatmypVVUVORtmz17thUUFGSlpKR42/bu3Wu53W7r0Ucf9Xl/p06drNjYWKuwsNDbNn36dCs4ONhKT08v8zumpKRYknw+o7aaPduyQkMtKzfX7pEAAAAUq668Vi0z2OHh4WrYsOE5+2VnZysxMVG33367wsPDve3jx49XeHi4z8z0ypUrVVBQoMmTJ/tcY/LkycrIyFBSUpK37Z133lGPHj0UGxvrbWvfvr0GDRrkc820tDSlpaVp4sSJCgoqvjVTpkxRUVGRVqxYUbEvXkutXCkNHSqFhdk9EgAAgOrnqG36du3apYKCAsXFxfm016lTRzExMdqxY4e3LTU1VeHh4erQoYNP3549e8qyLG9fy7K0c+dOde/evdTn9ezZU/v37/fWd+/YsUMul6vU50dFRenSSy/1+Xz4d+SItHmzb3kIAABAbeKogJ2ZmSmXy6WoqKhSr0VFRenw4cM+fZv5WS3nea+n77Fjx5SXl1fmNUv2zczM9Gk/2+fDv/ffNzXWv/613SMBAACwR4W36bMsS/n5+eXqGxISUqFrnzx5ssz31atXz/u6p29Z/Upe61zXrEjf7Ozs8n+ZWmrlSql/f7NjCAAAQG1U4YC9fv16DRgw4Jz9XC6X0tLS1K4Cp4qEhoZKkvLy8kq9durUKe/rnr5l9St5rXNdsyJ9S35+WaZOnarIyEiftrFjx2rs2LHnfG+gy86WEhOl2bPtHgkAAKjtli9fruXLl/u0ZWVlVctnVzhgd+jQQUuWLClXX3+lFufqb1mWt1SjpMzMTDVv3tyn77p16/z2k+Tt27hxY4WEhJR5zZJ9PePNzMxUdHR0qb69evU653eYP3++z2LK2uRf/5Ly86WEBLtHAgAAajt/E5zbt28vtdauKlQ4YDdr1kzjx4+virGoS5cucrvdSk5O1siRI73tp0+fVmpqqm655RZvW0xMjBYvXqy9e/f6LHTcvHmzXC6XYmJiJJmZ9K5duyo5ObnU523ZskVt2rTx7lgSExMjy7KUnJzssygyMzNTGRkZuuuuuyr9O9ckH34odewotW5t90gAAADs46hFjhERERo8eLCWLl3qc3Lja6+9ptzcXJ/DZhISEuR2u7Vw4UKfa7zwwguKjo5Wnz59vG0jR47Utm3btH37dm9benq6PvvsM59rdurUSR06dNCLL74oy7K87QsXLlRQUJBuvvnmSv2+NYllSR99JA0bZvdIAAAA7FXhGezzNXPmTLlcLu3Zs0eWZem1117Thg0bJEmPPvqot99TTz2lvn37qn///po4caIOHjyoefPmaejQoRoyZIi3X3R0tO677z7NnTtX+fn56tGjh9577z1t3LhRy5Yt8zktcsqUKVq0aJGGDx+uBx54QG63W/Pnz1dUVJTuv/9+n3HOmTNHCQkJGjJkiMaMGaNdu3ZpwYIFmjBhgtq3b1/Fdylw7dwpHT5MwAYAAKiWkxwty7JcLpcVFBRU6hEcHFyq78aNG61+/fpZYWFhVrNmzax77rnH52THkp5++mmrdevWVr169ayuXbtay5cv99vv0KFD1ujRo62GDRtaERERVkJCgrV//36/fVetWmXFxsZaoaGhVsuWLa0nnnjCKigoOOv3q+0nOc6aZVnh4ZZ16pTdIwEAAPCvuvKay7JK1ELgvHmK5lNSUmrlIsf+/aVGjaRVq+weCQAAgH/VldccVYONwPTzz9KmTZSHAAAASARsVILERKmwkIANAAAgEbBRCT78UOrUSWrVyu6RAAAA2I+AjQvC9nwAAAC+CNi4IF9+KWVmErABAAA8CNi4IGvXSvXrS/362T0SAAAAZyBg44J8+KE0aJAUEmL3SAAAAJyBgI3z9vPPUlIS5SEAAAAlEbBx3tauNdvzDR9u90gAAACcg4CN87Z6tRQbK7VoYfdIAAAAnIOAjfOSn29msG+4we6RAAAAOAsBG+dlwwYpK4uADQAAcCYCNs7L6tXSpZdKMTF2jwQAAMBZCNioMMsyAfuGGySXy+7RAAAAOAsBGxW2e7f03XeUhwAAAPhDwEaFrV4tNWggxcfbPRIAAADnIWCjwlavlq67jtMbAQAA/CFgo0IyM6WtWykPAQAAKAsBGxXywQdScDCnNwIAAJSFgI0KWb1a6tdPatzY7pEAAAA4EwEb5fbzz1JiIuUhAAAAZ0PARrm98opUWCiNG2f3SAAAAJyLgI1yKSqSFiyQRo2SLrnE7tEAAAA4l9vuASAwfPSRtH+/9Prrdo8EAADA2ZjBRrk8/7wUGyv17m33SAAAAJyNGWyc07590tq1pgbb5bJ7NAAAAM7GDDbOaeFCqUkT6ZZb7B4JAACA8xGwcVY5OdLLL0t33imFhto9GgAAAOcjYOOs3nhDys6WJk+2eyQAAACBgYCNs/rHP8zBMq1a2T0SAACAwEDARplycqTUVE5uBAAAqAgCNsq0a5dkWVK3bnaPBAAAIHAQsFGmHTukOnWkTp3sHgkAAEDgIGCjTDt2SJ07S3Xr2j0SAACAwEHARplSU6WYGLtHAQAAEFgI2PDr9GlTg039NQAAQMUQsOFXerqUl8cMNgAAQEURsOHXjh3m+aqr7B0HAABAoCFgw6/UVKlNGyky0u6RAAAABBYCNvzasYP6awAAgPNBwEYplsUOIgAAAOeLgI1SfvhBOn6cGWwAAIDzQcBGKamp5pkZbAAAgIojYKOUHTukiy6Smje3eyQAAACBh4CNUjz11y6X3SMBAAAIPARslMIOIgAAAOePgA0f//2vWeRI/TUAAMD5IWDDx5dfmmdmsAEAAM4PARs+duyQwsKktm3tHgkAAEBgImDDR2qqdOWVUnCw3SMBAAAITARs+EhKkuLi7B4FAABA4CJgw+vAAWn/fmnwYLtHAgAAELgI2PD69FMpKEiKj7d7JAAAAIGLgA2vxESpRw+pYUO7RwIAABC4CNiQJBUVmRlsykMAAAAuDAEbksz+10ePSkOG2D0SAACAwEbAhiRTHhIWJvXubfdIAAAAAhsBG5JMwO7fXwoJsXskAAAAgY2ADZ06Ja1fT/01AABAZSBgQ5s2mZBN/TUAAMCFI2BDiYnSxRdLXbrYPRIAAIDAR8CGPvlEGjTIHDIDAACAC0OkquWOHZNSUqi/BgAAqCwE7Fru888lyyJgAwAAVBYCdi23cqXUtq3UsqXdIwEAAKgZCNi12JYt0htvSH/4g90jAQAAqDkI2LVUQYF0111St27SlCl2jwYAAKDmcNs9ANhjwQLpyy/NLLabfwUAAACVhhnsWujQIWn6dGnyZKlHD7tHAwAAULMQsGuhqVOl8HDpqafsHgkAAEDNQ3FALfOvf0lvvy0tXSo1bGj3aAAAAGoeZrBrmblzpV69pHHj7B4JAABAzcQMdi3y3XfmWPRXXpFcLrtHAwAAUDMxg12LvPKK1KCBNGqU3SMBAACouQjYtURhofTyy9LYsWaBIwAAAKoGAbuW+OQTKSNDuuMOu0cCAABQsxGwa4mXXpK6dmXfawAAgKpGwK4F/vMfafVqM3vN4kYAAICqRcCuBV5/3QTr226zeyQAAAA1HwG7hrMsafFiacQIqUkTu0cDAABQ8xGwa7j166W0NBY3AgAAVBcCdg32n/9I48dLcXHSoEF2jwYAAKB2qJaAfeTIET388MMaOHCgIiIiFBQUpPXr1/vtGx8fr6CgoFKP4cOHl+prWZaeeeYZtWnTRqGhobrqqqv05ptv+r3u4cOHNXr0aDVq1EiRkZG68cYbdeDAAb99V69erbi4OIWGhqpVq1aaMWOGCgsLz/8G2CAvT7rpJvP83ntSEL9KAQAAVItqOSo9PT1dc+bMUdu2bXXllVcqKSmpzL4ul0stWrTQ008/LcuyvO3Nmzcv1XfatGmaPXu2Jk2apO7du2vVqlUaN26cgoKCNHr0aG+/3NxcxcfHKzs7W9OnT5fb7da8efMUHx+v1NRUNWrUyNt37dq1GjFihAYOHKjnn39eu3bt0syZM/XTTz9pwYIFlXRHqpZlSXfdJSUnS+vWSS1a2D0iAACAWsSqBjk5Odbx48cty7KsFStWWEFBQdYXX3zht298fLzVtWvXc17z0KFDVt26da177rnHp71///5Wy5YtraKiIm/b7NmzraCgICslJcXbtnfvXsvtdluPPvqoz/s7depkxcbGWoWFhd626dOnW8HBwVZ6enqZ40lJSbEk+XyGXebOtSzJsl57ze6RAAAAOEd15bVqKRwIDw9Xw4YNK/SewsJC5ebmlvn6ypUrVVBQoMmTJ/u0T548WRkZGT6z5O+884569Oih2NhYb1v79u01aNAgvfXWW962tLQ0paWlaeLEiQoqUVMxZcoUFRUVacWKFRX6DnbYuVN68EHpoYek22+3ezQAAAC1jyMrc7/++muFh4erQYMGioqK0uOPP66CggKfPqmpqQoPD1eHDh182nv27CnLsrRjxw5Jpk57586d6t69e6nP6dmzp/bv3+8N8jt27JDL5VJcXJxPv6ioKF166aXeazrZJ59I9epJf/6z3SMBAAConaqlBrsirrjiCg0cOFBdu3ZVbm6uVqxYoZkzZ+qbb77R8uXLvf0yMzPVrFmzUu+PioqSZBY1StKxY8eUl5fnbS+rb9u2bZWZmenTfmZfzzWdLCnJHIdep47dIwEAAKidKhywLctSfn5+ufqGhIRUeECLFi3y+fnWW2/VpEmT9NJLL2nq1Knq2bOnJOnkyZN+r1+vXj3v6yWfK6NvdnZ2hb9PdbIsE7ApDQEAALBPhQP2+vXrNWDAgHP2c7lcSktLU7t27c5rYCX93//9nxYtWqTExERvwA4NDVVeXl6pvqdOnfK+XvK5Mvp6Xj+bqVOnKjIy0qdt7NixGjt27Dnfe6EyMqTDh6Xevav8owAAABxt+fLlPtUPkpSVlVUtn13hgN2hQwctWbKkXH39lVqcjxb/22fu2LFjPtdet25dqb6eMg/Ptn6NGzdWSEiIt/1sfT3jzczMVHR0dKm+vXr1OudY58+f77OYsjp51nX+6le2fDwAAIBj+Jvg3L59e6m1dlWhwgG7WbNmGj9+fFWMpUz79++XJF100UXetpiYGC1evFh79+71Wei4efNmuVwuxcTESDIz6V27dlVycnKp627ZskVt2rRReHi495qWZSk5OdlnUWRmZqYyMjJ01113Vcn3qyxJSVLr1pKf0nQAAABUE0ftIpKdne23vnvmzJlyuVwaOnSoty0hIUFut1sLFy706fvCCy8oOjpaffr08baNHDlS27Zt0/bt271t6enp+uyzz3wOpOnUqZM6dOigF1980eeQm4ULFyooKEg333xzpXzPqrJ5M+UhAAAAdqu2XUQ8IXnPnj2yLEuvvfaaNmzYIEl69NFHJZlpe890/hVXXKGTJ0/q3XffVVJSkiZNmuSdlZak6Oho3XfffZo7d67y8/PVo0cPvffee9q4caOWLVsml8vl7TtlyhQtWrRIw4cP1wMPPCC326358+crKipK999/v88458yZo4SEBA0ZMkRjxozRrl27tGDBAk2YMEHt27evhjt1fvLypO3bpXHj7B4JAABALVelx9iU4HK5rKCgoFKP4OBgb58DBw5Yt9xyi9WmTRsrLCzMql+/vtWjRw9r0aJFZV736aeftlq3bm3Vq1fP6tq1q7V8+XK//Q4dOmSNHj3aatiwoRUREWElJCRY+/fv99t31apVVmxsrBUaGmq1bNnSeuKJJ6yCgoKzfj+7T3LctMmc3rhtmy0fDwAA4HjVlddcllWiFgLnzVM0n5KSYssix/nzpWnTpBMn2AMbAADAn+rKa46qwcb5S0qSuncnXAMAANiNgF1DJCWxPR8AAIATELBrgIwM8yBgAwAA2I+AXQNs3mye2aIPAADAfgTsGmDzZqlVK6mSDs4EAADABSBg1wDUXwMAADgHATvA5edLKSmUhwAAADgFATvAHThgTnG86iq7RwIAAACJgB3wcnLMc2SkveMAAACAQcAOcJ6AXb++veMAAACAQcAOcLm55pmADQAA4AwE7ADnmcEOD7d3HAAAADAI2AGOgA0AAOAsBOwAl5MjhYZKwcF2jwQAAAASATvg5eRQfw0AAOAkBOwAl5tLwAYAAHASAnaAy8mh/hoAAMBJCNgBjhIRAAAAZyFgBzgCNgAAgLMQsAMcARsAAMBZCNgBjkWOAAAAzkLADnAscgQAAHAWAnaAo0QEAADAWQjYAY6ADQAA4CwE7ABHwAYAAHAWAnYAsywWOQIAADgNATuA5edLBQUscgQAAHASAnYAy8kxz8xgAwAAOAcBO4ARsAEAAJyHgB3ACNgAAADOQ8AOYLm55pmADQAA4BwE7ADmmcFmkSMAAIBzELADGCUiAAAAzkPADmAEbAAAAOchYAewnBzJ5ZJCQ+0eCQAAADwI2AHghReknTtLt3tOcXS5qn9MAAAA8I+A7XAnTkh/+IP0xhulX8vJYYEjAACA0xCwHW79eqmwUDp+vPRrOTnUXwMAADgNAdvhEhPN87FjpV8jYAMAADgPAdvhPv3UPBOwAQAAAgMB28GOHJF275aaNvUfsD2LHAEAAOAcBGwH++wz8zxiRNkz2CxyBAAAcBYCtoMlJkpdukgdO1IiAgAAECgI2A5lWSZgDx4sNW5sykHy8nz7ELABAACch4DtUPv2SQcPSoMGmYAtld6qj4ANAADgPARsh0pMlIKDpWuuKQ7YZ5aJsMgRAADAeQjYDvXpp1KvXlKDBmUHbBY5AgAAOA8B24EKC80OIoMHm5/9BeyiImawAQAAnIiA7UCpqabeetAg83OjRua5ZMA+edIshCRgAwAAOAsB24ESE6WwMKl3b/Nz3bomSJcM2Dk55pmADQAA4CwEbAf64gvp6qtNsPZo3Ng3YOfmmmcCNgAAgLMQsB0oM1O6/HLftjMDtmcGm0WOAAAAzkLAdqCsLCkiwretrIDNDDYAAICzELAdKCtLioz0bSNgAwAABAYCtsNYlnTiBDPYAAAAgYqA7TAnT0oFBeeewWaRIwAAgDMRsB3mxAnzXJ4SEbfbd6cRAAAA2I+A7TBZWebZX4lIVpY55VEyAZvZawAAAOchYDuMJ2D7m8GWpJ9/Ns8EbAAAAGciYDuMp0TE3wy2VFwmQsAGAABwJgK2w5Q1g92okXn2BOzcXAI2AACAExGwHcYzg92ggW+7vxlsTnEEAABwHgK2w2RlSWFhUp06vu2UiAAAAAQGArbDnDhRujxEkkJDpZAQAjYAAIDTEbAdJiur9AJHSXK5fPfCJmADAAA4EwHbYbKy/M9gS74Bm0WOAAAAzkTAdpgTJ/zPYEulZ7BZ5AgAAOA8BGyHKe8MNiUiAAAAzkTAdpiyFjlKBGwAAIBAQMB2mLIWOUrFAbuwUDp1ioANAADgRARshylPiUhurvmZgA0AAOA8BGyHKc8iR89pjyxyBAAAcB4CtoMUFpra6rPNYBcVSZmZ5mdmsAEAAJyHgO0g2dnm+Wwz2JL0ww/mmYANAADgPARsB8nKMs9nm8GWCNgAAABORsB2EE9t9bkC9sGD5pmADQAA4DwEbAfxzGCXt0SERY4AAADOQ8B2kHOViDRoIAUHE7ABAACcjIDtIJ4SkbJmsF0uM4v9ww9SvXqS2119YwMAAED5ELAdJCvLzFCfbWa6cWPpxx+pvwYAAHAqAraDeA6ZcbnK7uOpwyZgAwAAOBMB20GyssouD/HwBGzqrwEAAJypWgL2Z599pjvuuEPt27dXeHi4Lr/8ck2YMEFHjhzx23/Tpk3q16+fwsPDFRUVpXvvvVe5ubml+lmWpWeeeUZt2rRRaGiorrrqKr355pt+r3n48GGNHj1ajRo1UmRkpG688UYdOHDAb9/Vq1crLi5OoaGhatWqlWbMmKHCwsLzvwHllJVV9gJHD2awAQAAnK1alsk99NBDOn78uEaNGqW2bdvq22+/1d/+9jetWbNGqampuvjii719U1NTNXjwYHXq1Enz589XRkaG5syZo3379mnNmjU+1502bZpmz56tSZMmqXv37lq1apXGjRunoKAgjR492tsvNzdX8fHxys7O1vTp0+V2uzVv3jzFx8crNTVVjRo18vZdu3atRowYoYEDB+r555/Xrl27NHPmTP30009asGBBld4nT4nI2RCwAQAAHM6qBhs2bCjVtn79esvlclmPPfaYT/uwYcOs6OhoKycnx9v20ksvWUFBQdYnn3zibTt06JBVt25d65577vF5f//+/a2WLVtaRUVF3rbZs2dbQUFBVkpKirdt7969ltvtth599FGf93fq1MmKjY21CgsLvW3Tp0+3goODrfT09DK/Y0pKiiXJ5zMq6oYbLOv668/e509/sizJshISzvtjAAAAaqXKyGvlUS0lIv369SvVdvXVV6tx48ZKS0vztmVnZysxMVG33367wksUGY8fP17h4eF66623vG0rV65UQUGBJk+e7HPdyZMnKyMjQ0lJSd62d955Rz169FBsbKy3rX379ho0aJDPNdPS0pSWlqaJEycqKKj41kyZMkVFRUVasWLFed6B8mEGGwAAIPDZtsgxNzdXOTk5atq0qbdt165dKigoUFxcnE/fOnXqKCYmRjt27PC2paamKjw8XB06dPDp27NnT1mW5e1rWZZ27typ7t27lxpDz549tX//fm99944dO+RyuUp9flRUlC699FKfz68KFanBZpEjAACAM9kWsOfPn6/Tp09rzJgx3rbMzEy5XC5FRUWV6h8VFaXDhw/79G3WrJnffpK8fY8dO6a8vLwyr1myb2Zmpk/72T6/KrDIEQAAIPBVeJGjZVnKz88vV9+QkBC/7evXr9eTTz6pW265Rddcc423/eTJk2W+r169et7XPX3L6nVRVbYAABzXSURBVFfyWue6ZkX6Zmdn+/0+laU8JSKe9ZgEbAAAAGeqcMBev369BgwYcM5+LpdLaWlpateunU/73r17ddNNN+nKK6/UokWLfF4LDQ2VJOXl5ZW63qlTp7yve/qW1a/ktc51zYr0Lfn5ZZk6daoiz5iGHjt2rMaOHXvW91kWM9gAAACVZfny5Vq+fLlPW1ZWVrV8doUDdocOHbRkyZJy9T2z1OLgwYO69tpr1ahRI61Zs8ZnIaOnv2VZ3lKNkjIzM9W8eXOfvuvWrfPbT5K3b+PGjRUSElLmNUv29Yw3MzNT0dHRpfr26tXrrN9XMqUvJRdTlldennT6NIscAQAAKoO/Cc7t27eXWmtXFSocsJs1a6bx48dX+IOOHTuma6+9VgUFBVq3bp3f+ukuXbrI7XYrOTlZI0eO9LafPn1aqampuuWWW7xtMTExWrx4sfbu3euz0HHz5s1yuVyKiYmRZGbSu3btquTk5FKft2XLFrVp08Yb9GNiYmRZlpKTk30WRWZmZiojI0N33XVXhb93eXl+oTrXDHbDhlK3blKXLlU2FAAAAFyAalnk+Msvv2jYsGHKzMzUhx9+qDZt2vjtFxERocGDB2vp0qU+Jze+9tprys3N9Tk8JiEhQW63WwsXLvS5xgsvvKDo6Gj16dPH2zZy5Eht27ZN27dv97alp6frs88+87lmp06d1KFDB7344ouyLMvbvnDhQgUFBenmm28+/5twDidOmOdzBezgYGn7dunqq6tsKAAAALgA1XKS47hx47Rt2zbdcccd2rNnj/bs2eN9rX79+kpISPD+/NRTT6lv377q37+/Jk6cqIMHD2revHkaOnSohgwZ4u0XHR2t++67T3PnzlV+fr569Oih9957Txs3btSyZcvkcrm8fadMmaJFixZp+PDheuCBB+R2uzV//nxFRUXp/vvv9xnrnDlzlJCQoCFDhmjMmDHatWuXFixYoAkTJqh9+/ZVdo88M9jnKhEBAACAw1XpMTb/c9lll1lBQUF+H61bty7Vf+PGjVa//9/enUdHWd1/HP88Q8hCEpYkLIkVCJtsoQEMWrESoQWpBhCRTQtHNgulUdoeY3FrgVClCByLBg0crQcXXFDpoW6IiscfWzDWLSBQE5MYCQqahRCy3N8fc2aSYSaQxMkMZN6vc3Imuc+dJ3c+HXm+vblzn6uvNu3atTNdu3Y1qampLnd2rO/BBx808fHxJjQ01CQkJJjnnnvOY7/CwkIzdepU07FjR9O+fXszceJEc/ToUY99X3vtNTNs2DATFhZmunfvbh544AFTXV19ztf4U+8MtGOH/Q6NDQwJAAAAP5Gv7uRoGVNvLQSazbFo/sCBA+f9kGNJiXT8uNS7d13bK69Ikyfb2+vdewcAAABe0pR67afw241mAtmaNVJysmsbS0QAAABaBwpsP8jPlwoKpO++q2srKZFCQ6XgYP+NCwAAAD8dBbYfFBfbHz/7rK6tMTeZAQAAwIWPAtsPjh+3P1JgAwAAtD4U2H7gmMH+9NO6tpIS1l8DAAC0BhTYfsAMNgAAQOtFge1jFRVSWZk0aJC9wHZsksgMNgAAQOtAge1jjtnr0aPtRXV+vv1nZrABAABaBwpsH3Osvx492v7oWCZSUkKBDQAA0BpQYPuYYwZ7+HApIqKuwP7xR5aIAAAAtAYU2D7mmMHu0kUaPLhuJxGWiAAAALQOFNg+dvy4faY6JERKSLDPYNfWSqWlzGADAAC0BhTYPlZcbJ+9luwz2Dk50g8/2H9mBhsAAODiR4HtY8ePS507278fPFiqrJQOHLD/zAw2AADAxY8C28fqz2AnJNgfP/zQ/sgMNgAAwMWPAtvH6s9gd+5sL7YpsAEAAFoPCmwfqz+DLdmXiezZY/+eJSIAAAAXPwpsH6s/gy3Zl4mUldm/ZwYbAADg4keB7UPl5dKpU+4z2JJkWVJ4uH/GBQAAAO+hwPYhx10cPRXY7dtLNv7XAAAAuOhR0vmQ4y6O9ZeIDBpkf2R5CAAAQOsQ5O8BBJL6t0l3iIyUevaUIiL8MiQAAAB4GQW2DzmWiMTEuLYnJNTdzREAAAAXNwpsHyouljp1ktq2dW1ftcr+4UcAAABc/CiwfejsLfoc+vf3/VgAAADQMviQow+dfZMZAAAAtD4U2D7U0Aw2AAAAWg8KbC8rLZWWLJF695a+/971GDPYAAAArR9rsL1s0iSpslKqqJD+7/+klJS6Y8xgAwAAtH7MYHvZyJHSkSP2mer9++vajWEGGwAAIBBQYHvZsmVSXJyUlORaYJeV2We2mcEGAABo3SiwW8jll9sLbGPsP3u6iyMAAABaHwrsFpKUZP+QY26u/WfHXRyZwQYAAGjdKLBbSFKS/dGxTIQZbAAAgMBAgd1CunSRunevK7AdM9jR0f4bEwAAAFoeBXYLqv9Bx+Jie3EdxMaIAAAArRoFdgtKSpIOHJBqatgDGwAAIFBQYLegpCT79nyHDrEHNgAAQKCgwG5Bw4fbH7Oy7DPYFNgAAACtHwV2C+rQQerXz74Ou7iYJSIAAACBgAK7hTk+6MgMNgAAQGCgwG5hSUnSxx8zgw0AABAo2DSuhSUlSZWV9u+ZwQYAAGj9mMFuYYmJUps29u+ZwQYAAGj9KLBbWLt20uDB9u+ZwQYAAGj9KLB9ICnJ/sgMNgAAQOtHge0DyclSp05SVJS/RwIAAICWRoHtAzNnSl99VbcWGwAAAK0XBbYPWJb9pjMAAABo/SiwAQAAAC+iwAYAAAC8iAIbAAAA8CIKbAAAAMCLKLABAAAAL6LABgAAALyIAhsAAADwIgpsAAAAwIsosAEAAAAvosAGAAAAvIgCGwAAAPAiCmwAAADAiyiwAQAAAC+iwAYAAAC8iAIbAAAA8CIKbAAAAMCLKLABAAAAL6LABgAAALyIAhsAAADwIgpsAAAAwIsosAEAAAAvosAGAAAAvIgCGwAAAPAiCmwAAADAiyiwAQAAAC+iwAYAAAC8iAIbAAAA8CIKbAAAAMCLKLABAAAAL6LABgAAALyIAhsAAADwIgpsAAAAwIsosAEAAAAv8kmBvXPnTs2dO1eXXXaZwsPD1bt3b82fP1/ffvutW9/k5GTZbDa3r9/85jdufY0xWrVqlXr16qWwsDD9/Oc/1/PPP+9xDN98842mTp2qTp06qUOHDpo0aZK++uorj323bdum4cOHKywsTD169NBf//pX1dTU/LQQAAAAEBB8UmCnpaXp/fff1+TJk/XPf/5TM2bM0AsvvKBhw4apuLjYpa9lWbr00kv1zDPPaPPmzc6vu+66y+28S5cu1d13361x48Zp/fr16tGjh2bOnKkXXnjBpV95ebmSk5P1wQcf6N5779WyZcuUnZ2t5ORknTx50qXv66+/rhtvvFFRUVFav369brzxRq1YsUKpqaneDwZ67rnn/D2EixK5NR2ZNQ+5NR2ZNQ+5NR2ZXcCMD3zwwQdubbt27TKWZZn77rvPpT05OdkkJCSc95yFhYUmODjYpKamurRfc801pnv37qa2ttbZ9tBDDxmbzWYOHDjgbDt48KAJCgoy99xzj8vzBw4caIYNG2Zqamqcbffee69p06aNOXToUIPjOXDggJHk8jtwfikpKf4ewkWJ3JqOzJqH3JqOzJqH3JqOzJrOV/WaT2awr776are2X/7yl4qKilJOTo7H59TU1Ki8vLzBc7766quqrq7WwoULXdoXLlyogoIC7d6929n28ssvKykpScOGDXO2XXbZZRozZozLbHdOTo5ycnK0YMEC2Wx10SxatEi1tbV66aWXzv9iAQAAEND89iHH8vJylZWVKSYmxu3Yl19+qfDwcEVGRio2Nlb333+/qqurXfp8/PHHCg8PV//+/V3aR4wYIWOMsrOzJdnXaX/yySe6/PLL3X7PiBEjdPToUWchn52dLcuyNHz4cJd+sbGx+tnPfuY8JwAAANCQIH/94rVr16qqqkrTp093ae/Tp49Gjx6thIQElZeX66WXXtKKFSt0+PBhl7VGRUVF6tq1q9t5Y2NjJdk/1ChJJ06cUGVlpbO9ob59+/ZVUVGRS/vZfR3nBAAAABrS5ALbGKMzZ840qm9ISIjH9l27dmnZsmWaNm2aRo0a5XIsMzPT5edbbrlFt99+uzZu3KglS5ZoxIgRkqSKigqP5w8NDXUer//ojb6lpaUNvNK65za05AWe/fjjj/roo4/8PYyLDrk1HZk1D7k1HZk1D7k1HZk1naNOc9RtLaXJBfauXbt07bXXnrefZVnKyclRv379XNoPHjyoyZMna8iQIW7FdEP+9Kc/KTMzUzt27HAW2GFhYaqsrHTre/r0aefx+o/e6Os47klubq4k6dZbb23Ua0Kds5fkoHHIrenIrHnIrenIrHnIrenIrHlyc3M1cuTIFjt/kwvs/v3766mnnmpU37OXWuTn52vs2LHq1KmTtm/frvDw8Ead59JLL5VkX+5R/9zvvfeeW1/HMo+4uDhJUlRUlEJCQpzt5+rrGG9RUZEuueQSt75XXHFFg2McN26cNm/erJ49e56zEAcAAIB/VFRUKDc3V+PGjWvR39PkArtr166aNWtWk3/RiRMnNHbsWFVXV+u9997zuH66IUePHpUkde7c2dmWmJioTZs26eDBgy4fdNyzZ48sy1JiYqIk+0x6QkKCsrKy3M67d+9e9erVy1noJyYmyhijrKwslw9FFhUVqaCgQL/73e8aHGNMTIxuueWWRr8mAAAA+F5Lzlw7+GQXkVOnTmn8+PEqKirSf/7zH/Xq1ctjv9LSUo/ru1esWCHLslz+38bEiRMVFBSkxx57zKXvhg0bdMkll+iqq65ytk2ZMkX79+93Wad06NAh7dy5U1OnTnW2DRw4UP3799cTTzwhY4yz/bHHHpPNZtNNN93U9BcPAACAgGKZ+pVkC5k0aZK2bdumuXPnKjk52eVYRESEJk6cKEl6//33NWPGDM2YMUN9+vRRRUWFtm7dqt27d+v22293K6bT0tK0evVqzZ8/X0lJSXrllVf0+uuv69lnn9W0adOc/crKyjR06FCVlpbqz3/+s4KCgrR27Vrndn7R0dHOvtu3b9fEiROVnJys6dOn69NPP9Wjjz6q+fPnKyMjo+VCAgAAQKvgkwI7Pj5eX3/9tcdjPXr00P/+9z9J9gXnd999t/bv369vv/1WNptNAwYM0IIFCzRv3jyPz3/ooYf0+OOPq6ioSH379tXSpUvdtv6T7FvxLVmyRG+99ZZqa2t17bXXas2aNR5n07dt26a//e1vysnJUefOnXXbbbfpvvvuU5s2bX5CCgAAAAgEPimwAQAAgEDhtzs5AgAAAK0RBfZPcObMGaWlpemSSy5Ru3btdOWVV2rHjh3+HtYFIysrS4sXL9bgwYMVERGhHj16aNq0aTp8+LBb34MHD+q6665TZGSkoqOjNWvWLH333Xd+GPWFJz09XTabTUOGDHE7Rm6uPvroI02YMEHR0dEKDw9XQkKC1q9f79KHzOocOXJE06dP16WXXqrw8HANGDBAy5cvd7sBQ6BmVl5ergceeEDjx49XdHS0bDabnn76aY99m5LRpk2bNHDgQIWFhalfv35u79GLWWMyM8boqaee0sSJE9W9e3dFREQoISFB6enpHu9DIbXuzKSmvdccqqurNXDgQNlsNq1Zs8Zjn9acW1MyM8YoIyNDQ4cOVbt27RQTE6MxY8bo008/devrtcwMmm369OkmODjYpKWlmczMTDNy5EjTtm1b8+GHH/p7aBeEKVOmmLi4OHPHHXeYTZs2mfT0dNOtWzcTERFhPv/8c2e/goICExMTY/r27WvWr19v/v73v5uoqCgzdOhQU1VV5cdX4H8FBQUmPDzcREZGmoSEBLdj5FbnzTffNCEhIeYXv/iFWbdundm4caP5y1/+YtLS0px9yKxOfn6+6dixo4mPjzcPPfSQyczMNHPmzDGWZZlJkyY5+wVyZrm5ucayLNOzZ08zevRoY7PZzL/+9S+3fk3JaMOGDcayLDN16lSzceNGM3v2bGNZllm1apWvXlaLakxmZWVlxrIsc9VVV5mVK1eajRs3mrlz55o2bdqY0aNHu52ztWdmTOPfa/U9/PDDJiIiwthsNvPwww+7HW/tuTUls9mzZ5vg4GAzb948s2nTJvPII4+Y2267zezYscOlnzczo8Bupr179xrLssyaNWucbadPnzZ9+vQxI0eO9OPILhy7d+92u7gcPnzYhIaGmt/+9rfOtoULF5rw8HBTUFDgbNuxY4exLMtkZmb6bLwXomnTpplf/epXJjk52a3AJrc6JSUlplu3bmbKlCnn7EdmddLT043NZjM5OTku7bNnzzY2m8388MMPxpjAzuzMmTPm2LFjxhhjsrKyjGVZHi/gjc2ooqLCxMTEmAkTJrg8/9ZbbzWRkZHOzC9mjcnszJkzZvfu3W7PXbZsmbHZbOadd95xtgVCZsY0/r3mcOzYMdOxY0ezYsUKY1mWW4EdCLk1NrMtW7YYy7LMa6+9ds7zeTszlog000svvaSgoCDNnz/f2RYSEqK5c+dq9+7dKiws9OPoLgxXXnmlgoJc72XUp08fDRo0SDk5Oc62rVu36oYbbnC5e+aYMWPUr18/vfDCCz4b74Vm165d2rp1q9atW+fxOLnVeeaZZ1RcXKz09HRJ9r33jYfPb5NZndLSUklSly5dXNq7desmm82m4OBgSYGdWdu2bd3y8aSxGb377rs6ceKEFi1a5PL83//+9yorK9P27du9N3g/aUxmbdu21ZVXXunWfuONN8oY43J9CITMpMa/1xzuvvtuDRgwoMEb3AVCbo3NbO3atbriiis0YcIEGWN06tQpj/28nRkFdjN9/PHH6tevnyIiIlzaR4wY4TwOz44dO6aYmBhJ9u0Ti4uLXe6c6TBixAhlZ2f7engXhNraWqWmpmr+/PkaNGiQ23Fyc/XOO++offv2ys/PV//+/RUREaH27dtr0aJFzjWdZOYqOTlZxhjNmTNH//3vf1VQUKAtW7Zow4YNuuOOOxQWFkZmjdCUjBzfDx8+3KXf8OHDZbPZAj7PoqIiSXJeHyQy82Tfvn16+umntW7dOlmW5bEPudmVlpZq3759SkpK0j333KMOHTooIiJCvXv31osvvujS19uZUWA3U1FRkWJjY93aY2NjZYzRN99844dRXfg2b96swsJC517ljn9QG8ryxIkTqqqq8ukYLwQZGRn6+uuvtXz5co/Hyc3V4cOHVVVVpYkTJ2r8+PHaunWr5s6dqw0bNmjOnDmSyOxs48aN0/Lly/X2229r6NCh6t69u2bOnKnU1FStXr1aEpk1RlMyKioqUps2bVwKSMk+ExcdHR3w141Vq1apQ4cOGj9+vLONzNz94Q9/0IwZM5wTep6Qm93Ro0dljNFzzz2nJ598UqtXr9azzz6rLl26aPr06Xrrrbecfb2dWdD5u8CTiooKhYSEuLWHhoY6j8PVwYMHtXjxYo0cOVKzZs2SVJfT+bJs27at7wbqZydOnNADDzyg+++/X1FRUR77kJursrIyVVRUaOHChVq7dq0k+x1kKysr9cQTT2jZsmVk5kHPnj01atQoTZkyRVFRUdq+fbvS09PVrVs3LVq0iMwaoSkZVVRUOJfeeOobyNeNlStXaufOncrIyFD79u2d7WTm6sknn9Tnn3+uV1555Zz9yM2urKxMkv26unfvXudfmlJSUhQfH68VK1Zo7NixkryfGQV2M4WFhXncTuj06dPO46hz7NgxXX/99erUqZNefPFF55+1HDmRZZ177rlH0dHRWrx4cYN9yM2V47WefRfXmTNn6vHHH9fu3bs1YMAASWTm8Pzzz2vBggU6cuSIc/Z10qRJqqmpUVpammbMmMH7rBGaklFYWJjOnDnj8TynT58O2Cy3bNmi++67T/PmzdOCBQtcjpFZndLSUi1dulR33XWX4uLiztmX3OwcrzM+Pt5lGVd4eLhSUlL0zDPPqLa2VjabzeuZsUSkmWJjY51/GqzP0Xa+N38gKSkp0XXXXaeSkhK98cYb6tatm/OY48LeUJZRUVEBNTt25MgRZWZmKjU1VYWFhcrLy1Nubq5Onz6tqqoq5eXl6eTJk+R2Fsd/b127dnVpd3wAhszcZWRkaNiwYW5LGyZMmKBTp04pOzubzBqhKRnFxsaqpqbGbX/sqqoqff/99wF53Xj77bc1e/ZspaSkKCMjw+04mdX5xz/+oaqqKk2dOlV5eXnKy8tTfn6+JPu/cXl5ec7lSORm19C1QbJfH6qqqlReXi7J+5lRYDdTYmKivvzyS+efHxz27Nkjy7KUmJjop5FdWCorK3XDDTfoyJEj2r59uy677DKX43FxcercubOysrLcnrtv376Ay7GwsFDGGKWmpio+Pl7x8fHq1auX9u7dq0OHDqlXr15avnw5uZ3F8aGUs3fvcayZ69KlC5md5dixY6qpqXFrr6qqkjFG1dXVZNYITckoMTFRxhi3vvv371dtbW3A5bl3715NnjxZI0aM0JYtW2SzuZckZFYnPz9fJ0+e1MCBA53Xh2uuuUaWZSk9PV29evVy7sBCbnaxsbHq1q2bx53dCgsLFRoaqsjISEktkFmTNvWDk2Mf7Pp7T1ZWVpq+ffuaq666yo8ju3DU1NSYCRMmmODgYPPGG2802O9ce8g+8cQTvhjqBeO7774zr732mtvX4MGDTc+ePc22bdvMZ599Zowht/qys7ONZVnm1ltvdWmfMWOGCQ4ONkVFRcYYMqsvJSXFhIaGmsOHD7u0T5o0yQQFBZHZWZq7D3b9jCoqKkx0dLTHfXYjIiLMyZMnW+4F+MG5Mvviiy9MTEyMGTJkyDn3Fw60zIxpOLfs7Gy3a0NmZqaxLMvMmTPHbNu2zZSUlBhjAi+3c73X7rzzTmOz2VxuKnP8+HHToUMHk5KS4mzzdmaWMR42i0WjTJs2Ta+++qruvPNO9enTR0899ZSysrK0c+dOjRw50t/D87s777xTjzzyiCZMmKCbb77Z7bhj/86CggINGzZMHTp00B133KHS0lKtXr1a3bt31759+wL+T9CSdO211+r777/XJ5984mwjN1fz5s3Tk08+qZtvvlmjRo3Su+++q5dffllLly517sZCZnU++OADjRkzRlFRUVq8eLGio6P173//W2+++abmz5+vDRs2SCKzRx99VD/88IMKCwu1YcMGTZ48WUOHDpUkpaamKjIyskkZZWRkaPHixbrppps0btw47dq1S5s3b9bKlSuVlpbmr5fpVefLzLIsDRw4UEVFRVq5cqXbn9579+7tsk92IGQmNe69dra8vDzFx8dr9erV+uMf/+hyLBBya0xmxcXFGjp0qMrLy7VkyRK1b99ejz/+uAoKCrRnzx4NHjzYeT6vZtakchwuKisrzV133WXi4uJMWFiYueKKK8zbb7/t72FdMJKTk43NZmvwq74vvvjCXHfddSYiIsJERUWZWbNmmeLiYj+N/MKTnJxshgwZ4tZObnWqq6vNsmXLTHx8vAkJCTH9+vUzjzzyiFs/Mquzf/9+c/3115u4uDgTEhJi+vfvbx588EFTU1Pj0i+QM+vZs2eD/4bl5eU5+zUlo40bN5oBAwaY0NBQ07dvX4/v04vZ+TLLzc0957Xhtttucztna8/MmMa/1+pzZOnpVunGtP7cGpvZV199ZW666SbTsWNHEx4ebn7961+bAwcOeDyntzJjBhsAAADwIj7kCAAAAHgRBTYAAADgRRTYAAAAgBdRYAMAAABeRIENAAAAeBEFNgAAAOBFFNgAAACAF1FgAwAAAF5EgQ0AAAB4EQU2AAAA4EUU2AAAAIAXUWADAAAAXvT/ofMU0r7/yDoAAAAASUVORK5CYII=", "text/plain": [ "PyPlot.Figure(PyObject )" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "1-element Array{Any,1}:\n", " PyObject " ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(elbo)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Moving forward with ReverseDiff\n", "\n", "- forward mode autodiff for $f: \\mathbb{R}^n \\rightarrow \\mathbb{R}$ is $\\mathcal{O}(n)$\n", "- reverse mode is only $\\mathcal{O}(1)$\n", "- reverse mode lies behind Theano, TensorFlow, Stan\n", "\n", "Once distributions are parametric, defining a new type that **just works** becomes easy\n", "- same(-ish) idea as ForwardDiff.jl\n", "- caveat: egregious memory usage" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false, "slideshow": { "slide_type": "skip" } }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "WARNING: using DataStructures.update! in module Main conflicts with an existing identifier.\n" ] }, { "data": { "text/plain": [ "backprop (generic function with 3 methods)" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "include(\"/Users/jmxp/code/rdiff/prototype6.jl\")" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "collapsed": false, "slideshow": { "slide_type": "-" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " 1.447133 seconds (8.04 M allocations: 679.577 MB, 6.46% gc time)\n", " 0.025645 seconds (155.60 k allocations: 5.186 MB)\n" ] } ], "source": [ "@time ∇L(stor, xx)\n", "∇L_r = grad(ELBO, length(xx))\n", "rstor = similar(xx)\n", "∇L_r(xx, rstor);\n", "@time ∇L_r(xx, rstor);" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# VinDsl needs your help!\n", "\n", "![](http://cdn.badassdigest.com/uploads/images/30118/the_pacifier02crop__index.jpg)\n", "\n", "Still finalizing API\n", "- docs\n", "- tests\n", "- better ideas!\n", "- Reverse mode!" ] } ], "metadata": { "celltoolbar": "Slideshow", "kernelspec": { "display_name": "Julia 0.4.6", "language": "julia", "name": "julia-0.4" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "0.4.6" } }, "nbformat": 4, "nbformat_minor": 0 }