From 9591a459e7c53e785ae9dde54e507803457fc14c Mon Sep 17 00:00:00 2001 From: jaimerz Date: Thu, 8 Jun 2023 16:44:02 +0100 Subject: [PATCH] first commit --- Lab.ipynb | 316 +++++++++++++++++++++++++++++++++++++ src/inference/Inference.jl | 4 +- 2 files changed, 318 insertions(+), 2 deletions(-) create mode 100644 Lab.ipynb diff --git a/Lab.ipynb b/Lab.ipynb new file mode 100644 index 000000000..3b79e7266 --- /dev/null +++ b/Lab.ipynb @@ -0,0 +1,316 @@ +{ + "cells": [ + { + "attachments": {}, + "cell_type": "markdown", + "id": "7f5caced", + "metadata": {}, + "source": [ + "# No Glue Code" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "1b21f376", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\u001b[32m\u001b[1m Activating\u001b[22m\u001b[39m project at `~/Cambdrige`\n" + ] + } + ], + "source": [ + "using Pkg\n", + "Pkg.activate(\"..\")" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "fe8d527d", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\u001b[36m\u001b[1m[ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mPrecompiling Turing [fce5fe82-541a-59a6-adf8-730c64b5f9a0]\n", + "WARNING: Method definition sample(Random.AbstractRNG, AbstractMCMC.AbstractModel, AbstractMCMC.AbstractSampler, AbstractMCMC.AbstractMCMCEnsemble, Integer, Integer) in module AbstractMCMC at /home/jaimerz/.julia/packages/AbstractMCMC/bE6VB/src/sample.jl:81 overwritten in module Inference at /home/jaimerz/Cambdrige/Turing.jl/src/inference/Inference.jl:210.\n", + " ** incremental compilation may be fatally broken for this module **\n", + "\n", + "WARNING: Method definition kwcall(Any, typeof(StatsBase.sample), Random.AbstractRNG, AbstractMCMC.AbstractModel, AbstractMCMC.AbstractSampler, AbstractMCMC.AbstractMCMCEnsemble, Integer, Integer) in module AbstractMCMC at /home/jaimerz/.julia/packages/AbstractMCMC/bE6VB/src/sample.jl:81 overwritten in module Inference at /home/jaimerz/Cambdrige/Turing.jl/src/inference/Inference.jl:210.\n", + " ** incremental compilation may be fatally broken for this module **\n", + "\n" + ] + } + ], + "source": [ + "using Random\n", + "using LinearAlgebra\n", + "using PyPlot\n", + "\n", + "#What we are tweaking\n", + "using Revise\n", + "using Turing\n", + "using MCMCChains" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "id": "4e582ff1", + "metadata": {}, + "source": [ + "## Model" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "10db58eb", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "funnel (generic function with 2 methods)" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Just a simple Neal Funnel\n", + "d = 21\n", + "@model function funnel()\n", + " θ ~ Normal(0, 3)\n", + " z ~ MvNormal(zeros(d-1), exp(θ)*I)\n", + " x ~ MvNormal(z, I)\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "6acf8fad", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "DynamicPPL.Model{typeof(funnel), (), (), (), Tuple{}, Tuple{}, DynamicPPL.ConditionContext{NamedTuple{(:x,), Tuple{Vector{Float64}}}, DynamicPPL.DefaultContext}}(funnel, NamedTuple(), NamedTuple(), ConditionContext((x = [1.2142074831535152, 1.23371919965455, -0.8480146960461767, 0.1600994648479841, 1.9180385508479283, -3.401523464506408, -0.0957684186471088, 0.6734622629464286, -3.2749467689509633, -1.6760091758453226, 1.9567202902549736, 0.1136169088905351, 0.11117896909388916, -0.5373922347882832, -0.12436857036298687, -1.2901071061088532, 1.702584517514787, -0.44460133117954226, 1.0818722439221686, 1.2208011493237483],), DynamicPPL.DefaultContext()))" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "Random.seed!(1)\n", + "(;x) = rand(funnel() | (θ=0,))\n", + "funnel_model = funnel() | (;x)" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "id": "f6a25810", + "metadata": {}, + "source": [ + "## Sampling" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "15adcfac", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "NUTS{Turing.Essential.ForwardDiffAD{0}, (), AdvancedHMC.DiagEuclideanMetric}(500, 0.95, 10, 1000.0, 0.0)" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "nadapts=500 \n", + "TAP=0.95\n", + "nuts = Turing.NUTS(nadapts, TAP)" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "79eb09bb", + "metadata": {}, + "outputs": [ + { + "name": "stderr", + "output_type": "stream", + "text": [ + "\u001b[36m\u001b[1m┌ \u001b[22m\u001b[39m\u001b[36m\u001b[1mInfo: \u001b[22m\u001b[39mFound initial step size\n", + "\u001b[36m\u001b[1m└ \u001b[22m\u001b[39m ϵ = 0.8\n", + "\u001b[32mSampling: 100%|█████████████████████████████████████████| Time: 0:00:01\u001b[39m\n" + ] + }, + { + "data": { + "text/plain": [ + "Chains MCMC chain (5000×33×1 Array{Float64, 3}):\n", + "\n", + "Iterations = 501:1:5500\n", + "Number of chains = 1\n", + "Samples per chain = 5000\n", + "Wall duration = 8.31 seconds\n", + "Compute duration = 8.31 seconds\n", + "parameters = θ, z[1], z[2], z[3], z[4], z[5], z[6], z[7], z[8], z[9], z[10], z[11], z[12], z[13], z[14], z[15], z[16], z[17], z[18], z[19], z[20]\n", + "internals = lp, n_steps, is_accept, acceptance_rate, log_density, hamiltonian_energy, hamiltonian_energy_error, max_hamiltonian_energy_error, tree_depth, numerical_error, step_size, nom_step_size\n", + "\n", + "Summary Statistics\n", + " \u001b[1m parameters \u001b[0m \u001b[1m mean \u001b[0m \u001b[1m std \u001b[0m \u001b[1m mcse \u001b[0m \u001b[1m ess_bulk \u001b[0m \u001b[1m ess_tail \u001b[0m \u001b[1m rhat \u001b[0m ⋯\n", + " \u001b[90m Symbol \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m ⋯\n", + "\n", + " θ -0.0223 0.7886 0.0522 413.6563 169.7823 1.0060 ⋯\n", + " z[1] 0.6063 0.7498 0.0108 5069.7387 3427.1499 1.0025 ⋯\n", + " z[2] 0.6079 0.7573 0.0115 4320.9620 3052.8164 1.0003 ⋯\n", + " z[3] -0.4327 0.7209 0.0098 5429.8996 3366.9500 1.0010 ⋯\n", + " z[4] 0.0797 0.7164 0.0083 7441.3963 3238.2467 1.0015 ⋯\n", + " z[5] 0.9583 0.7697 0.0160 2240.6673 3829.7072 1.0005 ⋯\n", + " z[6] -1.7165 0.9111 0.0312 807.5832 658.5686 1.0026 ⋯\n", + " z[7] -0.0382 0.7061 0.0074 9082.6991 3630.8661 0.9999 ⋯\n", + " z[8] 0.3428 0.7174 0.0075 9121.5605 3459.5329 1.0005 ⋯\n", + " z[9] -1.6422 0.8804 0.0312 723.1318 862.1615 1.0048 ⋯\n", + " z[10] -0.8389 0.7586 0.0155 2380.5553 3225.9573 1.0015 ⋯\n", + " z[11] 0.9841 0.7839 0.0204 1480.5433 3333.3508 1.0031 ⋯\n", + " z[12] 0.0610 0.7260 0.0073 9992.2345 3227.2357 1.0002 ⋯\n", + " z[13] 0.0502 0.7218 0.0074 9597.7597 3492.8413 0.9998 ⋯\n", + " z[14] -0.2643 0.6985 0.0076 8465.5278 3550.8716 1.0000 ⋯\n", + " z[15] -0.0588 0.6992 0.0069 10097.0979 3600.5524 1.0000 ⋯\n", + " z[16] -0.6426 0.7403 0.0111 4449.4546 3389.4120 1.0015 ⋯\n", + " z[17] 0.8502 0.7630 0.0170 2148.6973 3579.2884 1.0014 ⋯\n", + " z[18] -0.2243 0.7161 0.0075 9092.2632 3689.1139 1.0003 ⋯\n", + " z[19] 0.5554 0.7510 0.0100 6084.4692 2931.5269 1.0000 ⋯\n", + " z[20] 0.6082 0.7676 0.0120 4119.3689 3506.0779 1.0006 ⋯\n", + "\u001b[36m 1 column omitted\u001b[0m\n", + "\n", + "Quantiles\n", + " \u001b[1m parameters \u001b[0m \u001b[1m 2.5% \u001b[0m \u001b[1m 25.0% \u001b[0m \u001b[1m 50.0% \u001b[0m \u001b[1m 75.0% \u001b[0m \u001b[1m 97.5% \u001b[0m\n", + " \u001b[90m Symbol \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m \u001b[90m Float64 \u001b[0m\n", + "\n", + " θ -2.2065 -0.3868 0.0925 0.4924 1.1762\n", + " z[1] -0.7890 0.0781 0.5685 1.0967 2.1797\n", + " z[2] -0.7992 0.0925 0.5585 1.0925 2.1798\n", + " z[3] -1.9621 -0.8904 -0.4010 0.0559 0.9111\n", + " z[4] -1.3592 -0.3646 0.0753 0.5318 1.5406\n", + " z[5] -0.4356 0.4178 0.9148 1.4547 2.5775\n", + " z[6] -3.5450 -2.3266 -1.6935 -1.0793 -0.0270\n", + " z[7] -1.4690 -0.4906 -0.0301 0.4206 1.3792\n", + " z[8] -1.0701 -0.1323 0.3267 0.7906 1.8310\n", + " z[9] -3.4659 -2.2309 -1.6129 -1.0068 -0.0858\n", + " z[10] -2.3800 -1.3378 -0.7948 -0.2874 0.5115\n", + " z[11] -0.4251 0.4377 0.9476 1.4964 2.6353\n", + " z[12] -1.4033 -0.4094 0.0648 0.5427 1.5338\n", + " z[13] -1.3641 -0.4022 0.0432 0.4999 1.4786\n", + " z[14] -1.7348 -0.6979 -0.2440 0.1765 1.1310\n", + " z[15] -1.4742 -0.4882 -0.0684 0.3913 1.3590\n", + " z[16] -2.1922 -1.1132 -0.6176 -0.1330 0.7404\n", + " z[17] -0.5195 0.3007 0.8186 1.3553 2.4104\n", + " z[18] -1.6663 -0.6915 -0.2125 0.2407 1.1570\n", + " z[19] -0.8265 0.0431 0.5089 1.0252 2.1040\n", + " z[20] -0.8001 0.0751 0.5598 1.1130 2.1932\n" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "nuts_samples = sample(funnel_model, nuts, 5000)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "ee06e93d", + "metadata": {}, + "outputs": [], + "source": [ + "theta_nuts = Vector(nuts_samples[\"θ\"][:, 1])\n", + "x10_nuts =Vector(nuts_samples[\"z[10]\"][:, 1]);" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "9e196c65", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAArEAAAL3CAYAAACd2x1cAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABTEUlEQVR4nO3deXxU1f3/8XcCSYCQhGDCVkISAoqKsotQFVQQt1LAWr91I0qpG1YEtUFRxGoRt1pXXMGtiop8XXCBgqJ+RSqbVgsISGRHFkkCagLJ/P7gl3vPhLnDTLaZM3k9Hw8ejzN3zr33zCTih3M+93PifD6fTwAAAIBF4iM9AAAAACBcBLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLFAHcrJyVFcXJzi4uL0+uuve/YbNGiQ4uLiNGPGDL/j+fn5iouLU35+ftD7zJgxQ3FxccrJyZEkffTRR859w/lz++23O9esqKjQjBkzNHjwYLVq1UoJCQlq2bKljjzySA0dOlT33HOPCgsLq/fFhGn//v2aP3++brzxRvXp00ctWrRQQkKC2rRpo6FDh2rOnDme527cuFFPPPGE/vSnP6lXr15KSkpSXFyc/vjHP9ZoTLfffvsh31+TJk3UqlUrdevWTfn5+XrppZf0yy+/VOv65s+wTZs22rdvX8B+mzZtcvpFUuX3Yf4OSe7nONzvcE1U/v4H+9OiRYs6u3+0qfx7p77++wQipXGkBwA0FLfccouGDRumxo3r/j+7Nm3aaOTIkYccX7Fihb788ku1bt1aZ5555iHvd+/eXZK0b98+/eY3v9GHH34oSerZs6dOOeUUNWrUSN99953ef/99vf3222rWrJnGjBlTp59FkhYuXKjBgwdLOvjZTjrpJCUnJ+u///2v3n77bb399tv605/+pGnTph0SzM2aNUvXX399nY3N/C7Ly8tVVFSkVatW6bnnntNzzz2nsWPH6uGHH9b//M//VPse27dv1/3336/bbruttoYdk5KTk/W73/0u4HvNmjWr59EAqGsEsUA9aNasmb799ls9/fTTuvLKK+v8fl26dDlkVlc6OFv25Zdfer5v9vvwww/Vrl07vffeezr++OP93i8qKtKsWbPUtm3bWh55YPHx8TrvvPN03XXX6eSTT/Z7b+bMmbrooov05JNP6te//rUuvfRSv/dzc3N17bXXqmfPnurZs6deffVV3XXXXbU2Nq/vct26dbr99tv14osv6g9/+IN2796tq6++OuzrN23aVL/88ovuu+8+XXXVVcrMzKyFUcemjIyMoL/XAGIL6QRAPbjuuuskSXfccYd++umnCI/m8F555RVJ0qRJkw4JYCUpLS1Nl19+uc4666x6Gc9pp52m119//ZAAVpIuuOACZ6n6+eefP+T93/72t3rooYeUn5+v448/vl5mwiUpLy9PL7zwgm688UZJB38Hvvvuu7Cv065dO/3ud79TSUmJ7rzzztoeJgBYiyAWqAdnn322BgwYoK1bt+rvf/97pIdzWNu3b5cktWrVKsIjCU2PHj0kHcx/jTZ33XWX2rVrpwMHDlT7Z3/XXXepcePGmjZtmtavXx/WuQcOHNDTTz+tgQMHqmXLlkpKSlJubq6uuuoqz+/rjTfe0B//+Ed17dpV6enpatKkiXJzc3X55Zdr9erV1foMgfzrX//Sb37zG7Vu3VoJCQlKT09X586ddfHFF+vjjz+utfuYCgsL/fLHA/HKKTWPf/jhhzrjjDOUnp6upk2bqmfPngH/ESVJAwcOVFxcnD766COtWLFCI0aMUEZGhpKSknTMMcfo/vvvl8/n8xzP/PnzNWLECLVt21aJiYlq1aqVhg8frkWLFlXnKwBiBkEsUE+mTp0qSbrnnnu0a9euCI8muA4dOkiSpk2bptLS0giP5vDWrFkjSfWW3hCOhIQEXXDBBZKkefPmVesanTt31ujRo1VWVqaJEyeGfF5JSYkGDx6s0aNHa+nSpTr++OM1dOhQJSUladq0aerRo4eWL19+yHm///3v9fLLL6tp06Y67bTTNGTIEMXHx2v69Onq1auXPvvss2p9DtNzzz2nM844Q3PmzFFubq7OO+88nXLKKUpNTdUrr7yiN954o8b3qCvPPvusTj/9dO3evVtnnnmmunfvruXLl2vkyJF68MEHPc/74IMP1LdvX61atUqDBw9Wv3799O233+qGG27wzNu+4YYbNGjQIL355pvq0KGDhg0bpo4dO+rNN9/UySefrOnTp9fRpwQs4ANQZ7Kzs32SfJ988onP5/P5RowY4ZPku/766/36nX766T5JvunTp/sdHzlypE+Sb+TIkUHvM336dJ8kX3Z2dtB+kyZN8knyDRgwIGi/v//97z5JPkm+1q1b+0aPHu175plnfMuWLfMdOHAg6Ln1bevWrb60tDSfJN9DDz102P6V38GoUaNqdN9Qv0ufz+d78cUXne9z//79IV3/ww8/9Eny5eXl+Xy+g58zOTnZFxcX51u+fLnTb+PGjc61q7rwwgt9knznnnuub/v27X7vVf6MO3fufMjP9JVXXvHt3bvX71hFRYXv0Ucf9UnyHXvssb6Kigq/9yu/j0mTJoX0+XJzc/3+2zBt377dt2zZspCu4/OF/vvv8/l869evP2zfyv9u169fH/B4QkKC7+233w44hrS0NN9PP/3k996AAQOcn9G0adP83ps/f74vLi7O16hRI9/GjRv93nvyySd9knydOnXyffnll37vLVy40JeSkuJLTEz0ffvttyGNH4g1zMQC9ehvf/ubGjdurMcee0zff/99pIfjaezYsbrrrruUnJys7du366mnntKoUaPUs2dPpaena+TIkbW6rFxdBw4c0MUXX6yioiIdd9xxuuKKKyI9pIAyMjKc9u7du6t1jTZt2uj666+Xz+dTQUHBYfuvXLlSL7/8stq1a6d//vOfh6SGjB07VmeffbbWrFmj9957z++9Cy64QMnJyX7H4uLidPXVV6tfv3765ptvtHLlymp9jkrbt29XWlqaTjrppEPea9WqlZMiEo7vv//es8TWRx99VKPxmq699lqde+65fsfy8/PVpUsXFRUVacmSJQHPGzFixCG/o5Uz3eXl5U41EOlgibvKcmWvvPLKIbnpp5xyim699VaVlZXpiSeeqIVPBdiHIBaoR0cddZQuv/xylZaW6tZbb430cIK6+eabtWnTJs2YMUOXXXaZunXrpkaNGqmkpETPP/+8evTooXfffTeiY7zyyis1f/58HXHEEXr99deVmJgY0fF4qaiocNo1qed64403KiMjQx988IFfwBPIu+++K5/Pp7POOkspKSkB+wwcOFCSAqYHrF27Vo888ojGjh2rUaNGKT8/X/n5+U6+dE3/EXPCCSeoqKhIl156qZYuXer3HVVXcnKyRo4cGfBPmzZtanz9Sr/5zW8CHj/66KMlSZs3b67xecuXL9eWLVuUl5enXr16BTwv2M8PaAgosQXUs8qySy+99JJuuOGGgE//V6oMeHxBHvow36/tgvctWrRwggBJ+vHHHzV79mxNnDhRW7du1ciRI/X9998ftgbnzp07dcMNNxxyvEuXLiHNKgZy3XXX6ZlnnlF6errmzZunI488slrXCeTpp5/Wp59+esjxgoICdenSJezr7dy5U9LBn096enq1x5WamqqJEydq7Nix+stf/qLFixd79q2shPDMM8/omWeeCXrdHTt2OO3y8nKNGTNGTzzxRNDfu+Li4jBH7++xxx7TueeeqxdeeEEvvPCCUlJS1KdPH5122mm65JJLnLzscNRXiS2vsaWmpkqS5wYX4ZxX+fNbt27dYf+7Nn9+QENCEAvUs7Zt2+q6667TlClTNGHChKC7TVUu6Xrt1lRp7969kqTmzZvX3kADSE9P1+WXX64ePXqoZ8+e2rlzp/7v//7P2Ygg2Piee+65Q44PGDCgWkHs+PHj9dBDD6lFixaaO3dutZaeg/n0008DjrdyyThcy5Ytk3QwaK9pia+rrrpKDz74oL744gu9/vrr6tevX8B+lTOb3bt3V7du3YJes2/fvk77H//4h6ZNm6Y2bdrogQceUP/+/dW6dWs1adJEknThhRfq5ZdfPuw/rA7n6KOP1urVqzV37lwtWLBAn332mT755BMtWLBAd9xxh5555hldfPHFNbpHdR1uVjg+vnqLmOGcVzmGNm3aaMiQIUH7mukqQENCEAtEwF/+8hc9+eSTevfdd4OWEqqcuVm7dm3Q61U+nV+d2avq6NGjhzIyMrRz505nljGYnJycGgc9lW666SY98MADSktL09y5c9W7d+9aua5pxowZtTajt3//fr366quSpDPOOKPG10tMTNRf//pXXXLJJbrllls0d+7cgP2ysrIkSb/+9a/1yCOPhHz9yrE+8cQTGjp06CHvV/6u1YbGjRvr7LPP1tlnny3p4OzuAw88oMmTJ+uKK67Q8OHDD8nNranKlJOSkpKA7+/fv19bt26t1XtWR+XP74gjjmADB8ADObFABKSlpenmm2+WdDAo83LaaadJkr766ivPQHb//v166623/PrX1OECzj179jjLye3bt6+Ve4aioKBA9957r9LS0jRv3jz16dOn3u5dXbfccou2bNmihISEWtv+9qKLLlK3bt20Zs0aPfXUUwH7VG5E8dZbb3kubwdS+eBZdnb2Ie998803WrFiRfgDDlFqaqpuv/12tWjRQj/99JO+/fbbWr9HZmamEhMTtXv3bv3www+HvP/BBx/owIEDtX7fcPXp00cZGRn673//q2+++SbSwwGiEkEsECHXXHONOnTooMWLF3sWLe/du7dOP/10+Xw+XXzxxYfMEP3888+66qqr9P333yszM1OjRo2qlbGdcMIJeuyxxwI+Sb9t2zaNHDlSZWVlys7O9lzOrm0TJ07U1KlT1aJFCysC2O+++06XXnqp7r33XknSI488EjAwrI64uDhNmTJFkjzrkvbo0UPnnXeeNm7cqBEjRhxSuF86mKby0ksvOQ9rSe5DRo8++qjfsvrWrVt16aWX1kqA99NPP+mBBx4ImMv5ySefaM+ePWrUqFGd/AMpISFBp5xyiqSDv1PmZ/zyyy81ZsyYWr9ndSQkJGjSpEny+XwaPnx4wBzt8vJyLViwQJ9//nkERghEHukEQIQkJSXpjjvuUH5+ftCtaF988UWdccYZWrx4sTp27Kj+/furbdu2Kioq0qJFi7Rr1y61bNlSs2bNUosWLWplbGvWrNE111yjP//5zzruuOOUl5enxo0ba/PmzVq8eLH279+vli1b6pVXXqmXbVzfeust3XXXXZKkTp066dFHHw3YLyMjQ/fdd5/fsa1bt2r48OHO602bNjnXPPHEE53jjz32mHr27Bn22FatWuVse1tRUaGioiKtWrVKa9askc/nU2Zmph555BH9/ve/D/vawZx11lkaOHBg0NJR06dP1549e/Tee+/pqKOOUrdu3ZSbmyufz6fCwkJ9+eWXKisr08qVK9W6dWtJB6tSvP/++3rqqaf04YcfqmfPniouLtbChQvVsWNHDR8+XLNnz67R2MvKyjR+/HjdeOONOu6449S5c2clJCSosLDQCchuueUWZWZm1ug+Xu688059/PHHeuqpp7Rw4UIdf/zx2rx5s5YsWaILL7xQH330UVSUwBszZow2bNige++9VyeffLKOPfZYderUSU2bNtW2bdu0YsUK7dmzR48//rjf7zLQUBDEAhF0ySWX6P7779d//vMfzz5t2rTR4sWL9cwzz2jWrFn68ssv9fHHH6tZs2bKy8vTn/70J/35z3+u1RJCn376qf71r39pwYIFWrNmjebPn6+9e/cqNTVVffr00ZAhQ3T11VfX2wMl5ozwkiVLPOtwZmdnHxLElpaWBnyKf8eOHX4zgdV92n779u3OQ2CJiYlKTU1Vu3btdMkll2jIkCEaMWKE81BUbZs6darfQ1lVpaSkaO7cuZo5c6ZefPFFLV26VCtWrFBqaqratm2riy66SEOHDlVeXp5zTt++fbVkyRJNnDhRX3zxhd566y1lZWXp2muv1cSJE3XttdfWeNzNmzfXtGnTtHDhQi1fvlzz5s1TWVmZ2rVrpxEjRujqq6+utdSYQPr27auFCxdq0qRJ+vzzz7Vx40YdeeSR+sc//qErr7xSubm5dXbvcN1zzz0aNmyYHnvsMX366ad6//33lZiYqLZt22rgwIE699xzNWLEiEgPE4iIOF9tPW0BAAAA1BNyYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYp3GkB1CfKioqtGXLFqWkpCguLi7SwwGAGvP5fCopKVG7du0UH8+8BICGo0EFsVu2bFFWVlakhwEAtW7jxo1q3759pIcBAPWmQQWxKSkpkg7+ZZ+amhrh0QBAzRUXFysrK8v5+w0AGooGFcRWphCkpqYSxAKIKaRIAWhoSKACAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFincaQHACCycgrmOO3Cu8+J4EgAAAgdM7EAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArNM40gMAYI+cgjlOu/DucyI4EgBAQ8dMLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrsO0sAAfbygIAbMFMLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrUGILiEG1XSrLvB4AANGAmVgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHWoEwvECGq5AgAaEmZiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2qEwANhFm9oPDucyI4EgAAao6ZWAAAAFiHIBYAAADWIZ0AaIDqcmME0hYAAPWBIBZAtRCsAgAiiXQCAAAAWIcgFgAAANYhiAUAAIB1yIkFEFBdPvwFAEBNMRMLAAAA6xDEAgAAwDoEsQAAALAOQSwAAACsQxALAAAA6xDEAgAAwDoEsQAAALAOQSwAAACsQxALAAAA6xDEAgAAwDpsOwvEOLaPBQDEImZiAQAAYB2CWAAAAFiHdALAYqQKAAAaKmZiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWoU4sgBqjXi0AoL4xEwsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrUJ0AsAyVAAAAYCYWAAAAFiKIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1qE6AYA6Z1ZUKLz7nAiOBAAQK5iJBQAAgHUIYgEAAGAd0gkAS7DJAQAALmZiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAddjsAECdYYMGAEBdYSYWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1mkc6QEAkHIK5jjtwrvPieBIAACwA0EsEMXM4BYAALhIJwAAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANZpHOkBAGi4cgrmOO3Cu8+J4EgAALZhJhYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB2CWAAAAFiHIBYAAADWIYgFAACAddh2Fogy5lassSjWPx8AoH4wEwsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrUJ0AQFQwqxYU3n1OBEcCALABM7EAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOtQJxaoR9RCBQCgdjATCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArNM40gMAYklOwRynXXj3OREcCQAAsY2ZWAAAAFiHIBYAAADWIYgFAACAdQhiAQAAYB0e7AIixHwIDP54QA4AcDjMxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOuw2QFQD9jYAACA2sVMLAAAAKxDEAsAAADrEMQCAADAOgSxAKJaTsEccooBAIcgiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcgFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWKdxpAcA2C6nYE6khwAAQIPDTCwAAACsQxALAAAA6xDEAgAAwDrkxAJ1hFxZAADqDjOxAAAAsA5BLAAAAKxDEAsAAADrEMQCAADAOjzYBVQDD20BABBZzMQCAADAOszEAmFgBhYAgOjATCwAAACsQxALAAAA6xDEAgAAwDoEsQAAALAOQSwAAACsQxALAAAA6xDEAgAAwDoEsQAAALAOQSwAAACsQxALAAAA6xDEAgAAwDoEsQAAALAOQSwAAACsQxALAAAA6zSO9AAAIBQ5BXOcduHd50RwJACAaMBMLAAAAKxDEAsAAADrEMQCAADAOuTEAodh5mICAIDoQBALwGo88AUADRPpBAAAALAOQSwAAACsQzoBAOuQpwwAYCYWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUosQUgZrB7FwA0HMzEAgAAwDoEsQAAALAO6QRAAOwIBQBAdGMmFgAAANYhiAUAAIB1CGIBAABgHYJYAAAAWIcHu4D/j4e5AACwBzOxAAAAsA5BLICYlFMwh9l1AIhhBLEAAACwDjmxaHDM2bnCu8+J4EgAAEB1MRMLAAAA6xDEAgAAwDqkE6DB4CEfAABiB0EsYhqBKwAAsYl0AgAAAFiHIBYAAADWIZ0AQINBeTUAiB3MxAIAAMA6zMSiQePBLwAA7MRMLAAAAKxDEAsAAADrkE4A61SmAHg9mEOKAEz8PgBAbCKIBdAgUakAAOxGOgEAAACsY81M7JQpU/TGG29o1apVatq0qfr376+pU6fqqKOOivTQEAVYMkZNMCsLAPaxJohduHChrrnmGvXp00cHDhzQzTffrDPOOEP//e9/lZycHOnhIQIIXFHXDpd/DQCIHGuC2Pfff9/v9YwZM9SqVSstXbpUp5xySoRGBVNdzmYRsAIAAJM1QWxVRUVFkqSWLVtGeCQNU3VnqAIFo+Y1CFYRaYF+B2vj99Lr95xZXgCoHiuD2IqKCo0dO1a//vWv1bVrV89+paWlKi0tdV7/+OOPkqRNmzYpNTW1zscZyw4U75R08Luseqzq8UDnmbyuAcSScP9bCdWePXskSbt3767RdQAgWvh8PpWUlKhdu3aKj/euQRDn8/l89TiuWnHVVVfpvffe06effqr27dt79rv99ts1efLkehwZAAAAasPGjRuDxnnWBbFjxozRm2++qY8//li5ublB+waaic3JydFJOluNlVDXQwWAOvez9mmx/qWNGzeywgQgJhQXFysrK0t79uxRWlqaZz9r0gl8Pp+uvfZazZ49Wx999NFhA1hJSkpKUlJS0iHHGytBjeMIYgHYL8GXKElKTU0liAUQU+Li4oK+b00Qe8011+if//yn3nzzTaWkpGjbtm2SpLS0NDVt2jTCowMAAEB9smbHrscff1xFRUUaOHCg2rZt6/yZOXNmpIcGAACAembNTKxlqbsAAACoQ9bMxAIAAACVCGIBAABgHYJYAAAAWMeanFgg1jVu3cppH9j+QwRHEjl8BwCAUDETCwAAAOsQxAIAAMA6pBMA9cxrydxr+bwhLbHH+ucDANQeZmIBAABgHYJYAAAAWId0AjQY9bEsH8o9wr232d+8fnWu25BSEwAAsY2ZWAAAAFiHIBYAAADWIZ0ADUa0L597LfXXNIXAqz+pBQAAmzETCwAAAOsQxAIAAMA6pBMANVTXy/JlXdo77fiFy2rtuvW9uQLpCwCA2sRMLADEgK6TPoj0EACgXhHEAgAAwDqkEyDm1GTZujrn1tbSuNd1EldtcvvUyp2qN45IXrdR16P8Xpd/vbqmwwEAWI6ZWAAAAFiHIBYAAADWIZ0AMac+UgjqWn2MyVyij/bl+WgfHwCg/hHEAkCMyCmY4/e68O5zIjQSAKh7pBMAAADAOgSxAAAAsA7pBIgJZg6pKZR8UrNPtOeien1OKbQduHyZ6U47bsePYd07VF7fYbjHAQAIhplYAAAAWIeZWACIUeaDXjzkBSDWEMQiJtRkGbomy9/VGZO5jG+mFpjHvc6tetw8v7FHP79zjHbFgJ5OO74ay/teaRGh/CxIIQAA1BTpBAAAALAOM7EAEKNIIQAQywhiUeeiZenYaxyhpAqEWhUglFQDr+V981yvigJVqxmEW2Hg5+F9nXbzNXvcN4zUAJ95gjE+M32g6ljCTb2oaXWCmlSjAADEBtIJAAAAYB2CWACIQaQSAIh1pBOgzkVLCkG4faozbnPpv7HHcTMFwDxuLuP7LdWHMFbJf7xVl/4rpXy2Pqxzva5TdSxlXdo77USP64Za1SHQvaumTdTkugCA2MBMLAAAAKzDTCwAxCBzo4NwkIYAwBYEsWgwavJEfLD+fsvZHhsWVJjL7UYfM23Ab+OCEKsTmBsWlPbPddpeaQNeaQ3mvfd2buG0k3bvd9oJu/b5Xcv8fOZmCTLG3iiEqgde33/Vz2oyr3sgSL9AoqVaBgCgZkgnAAAAgHWYiQUAOIKlIZBqACCaEMSiQfJKLTCfso9fuMxpe1UdqMrst+3CTk671bJfAvZpFELVAq9NBSTpgDHGpuY4PJ7s339EstNOMPp7pRCYqlYIMFMZEldtCjh2r80YzHQAr7SBoFUHwtzkwUQKAQDEBquC2I8//lj33nuvli5dqq1bt2r27NkaNmxYpIcFAA0CD4sBiCZW5cTu27dP3bp106OPPhrpoQAAACCCrJqJPeuss3TWWWdFehiIcl7L0Obyt5kqYC4vm4X6zafszSfrzeXyqszl+rTC8oD32+8xDq/i/n4bDlRZRl/34IlO+8inA6cNFJ3Q0mln/Hu30zZTCMxqBqGkWkj+1Qq8Ugi8NlEIliIRSLAUAK+fKwAgtlkVxIartLRUpaWlzuvi4uIIjgYAGqZw0xBIPwAQCqvSCcI1ZcoUpaWlOX+ysrIiPSQAAADUgpieiZ0wYYLGjRvnvC4uLiaQraFoLRTvNS7zeChLzV6fqVGQFILNg49w2mYVAr8leqO/+SS/mbJQZqQAmBsi+FUOMPpI0lFTvnPaJcZmB6bStLiAx702RPh5eF+n/Uu6++/c9GcX+XcMMz3Aq1KBOe6msxcH7BOM1881Wn9XAQC1I6aD2KSkJCUlJUV6GACAMHilH5BmAMAU0+kEAAAAiE1WzcTu3btXa9eudV6vX79eK1asUMuWLdWhQ4cIjiz2eC3FhrosG+5SbijpAMGuU1vLxeb9QnniXpLa7wi8fB7KmMx7mE/77zgnz2k3+bEiYJ+q55ubFBTluSsQv5q3y2mbn8Ncxm++Zo/TNtMMUox7mSkRVZnfm396gEd6hpGKYKYQhPrzDqUfKQQAENusCmKXLFmiU0891Xldme86cuRIzZgxI0KjAgDUh+putuCF9ATAblYFsQMHDpTP5zt8RwAAAMQ0q4JY1B+vpdi6Wt4Pd0m46hPwXv3M5XazwoD5NL0Xc8OAUr+NAfz7eVUFqNovEK/lfbMSwI+X93NPMMYRjHn+XrPagLHxQeacdU7b3Mhg0/lutYWOj7npO2bVgqqaG20zHcFnpA1sMqo4tP+ne12vzRy8Ujsk6YDHz6+maSkAAHsQxAIAGqTqpieQhgBEB6oTAAAAwDrMxCIsNizF+i0dh1CEPxTmk/9Vl7a9Ng0wmakC5rWqVhuoZC6xp61zt05el+//786UrxLd6xa5+eJJA3o6bb/KAx4bDpS2THDamSvKnbaZZlD1c5rvmakX8cZ3XmH0MaskeDF/v8w0iqobLfilIJg8Ph+pBaguZl2B6MVMLAAAAKxDEAsAAADrkE4AKwVbBjaXuePD3KjBa/MBs4qAuRGBJKWtS1Ygias2BTz+Q88mTvtX89x7/JJu/JvSqEJQlNPIaad85X+ttEJ36X/TULe982T3WkcXBK6A4MVMdzDTDIqqfO7StDin3WrZL067wkhlML/Dco/UDvO7NdMEMv692+1TNRXEI20gFKGkFtRmygHpC9GLVAHAbszEAgAAwDoEsQAAALAO6QSIOjVd4vVLAwihIoHPYxMEc1m8yNgkoCozbeC7qzs57dxd6YG6+z2lbz7V3+THioD9zaoDJdlxfu+ZKQjxRe57yd+7xwsfc7+DI14MeAu/lIW0wF0OYaYy+F0rL8l9YbSbdHY3SzBTFkxeKRhVef3MvJhpCqH0r81lf1IIoldtb2Nb30iHQEPHTCwAAACsQxALAAAA6xDEAgAAwDrkxCLqeOUQeuXHVn3PK+fRa8euOI9yTWZubYbRNvNYJf/8THO3K7OfmXP6S7qbX2vuRPWLsUOVySxlFUxCkftvUrPk1brj3ZJeZg6teV0z79bss/t493jLKuW9zJ3Adho5w5lz1jlts9yZWa7Lb9zGd2v+HL1KdUlVdggzjpvnxC9c5rS9ynt5/a7VR1ksSm8BQM0wEwsAAADrMBMLAICForm6ApUTUB8IYhF1QllmDbb86pU24HW+2d8sxbTX2DXLLAlVdWnbTEdI2RF4dyxzeX/98ESnnWYsf+882b1H3gy33Ja5DL+ju5uWIEnH/ukbp/3J58c47U2nuSkEGZ+4KQFmekBCkZmm4Lb3ZVcYfdzFmiY/+pfUMpfoS87v7/YzPreZclDasoXRTgh4PGWH+7MokzczVeDn4W7prqazFwfsH8rvhKk+lvercw9SEADARRALAABqVaRniZkJbhjIiQUAAIB1mIlF1KnOLl3mU+kKcdenQEJJDUjY5XeKdpyTF/BaGf/e7bS//aNbwcDcTcvc3SrjE/fcdfluakH7t7zHu+kv7g5hLfPclICdJ7uL8fFd3VSGIye51/1+sruknzQnxWkfcexWp93kejctYWeVXct2THZTCMyqDGYlhi0D3M991JTvvD/IYVStIGGmfZgpC+YObWbVCBnnB6tyEahPNInWcQFAJBDEAgAA65AyANIJAAAAYB1mYhHVQilGL0kHjKfVD4Rw3cYey85eGyWkfOa2zQL+ktTkxwoFYhbkb/mVWQnArRBQkn34jQy2DDD7+N/LawOB+CL3P+3GX6U67VXj3TSDpK/dKgmlacY13/2V+2Kw2zQrLPz/KzgtM4WgrEXAIfmlXZgbLaQVuqkIZhUIvwoGRjqH5J9C4MVrEwsz9cT8vTF/J6gCAADRjyAWAABYJ1IVEEhjiB6kEwAAAMA6zMQiqplPoZtL/aEu8Zrnx3k8oV5hpAckeNzPZBbal6SmRvvHy/s57SS3OIFfpYLNg49w2u0XuEv0/hsUuEvpaetKA/aRpC0D3PSCijR3Kd5MJ0gqctMXyja4KQTmvX/o6V635Hg35aD9W26agFlJQfLfOOHIp93Pt2q8W+mgy/0lCsSsdGCmEJjMjR3M6geS/8/SrByR4vEzNn8PEo3qFWbqCWkDAGAXglgAAIAQhZrGQNpB3SOdAAAAANZhJhZRzWtJvzrnexW291te9iiEb1YwaGQW0Zd/FQJz6d+8rnm+ubxvpgekfO8eNysemE/pV11WNysX5LTf4bQ3FLUNeK2dJ7vnr8t3/w2b8pV7TTMVYdNQd8G9aeo+v3snf+pWPdjbuYXTzvgkPuDxX9Ld4ztPdlMISrLd7yD3NXcniY6PBf65SJLPY7MDk5lC4Heu8bPw+gsw2IYIpB0AQHQgiAUAAKhltVk9gdSEwEgnAAAAgHWYiUVUC7XofCj9vJ5WPxBCyoGZMlB1gwFzOdtcPjcL9Kd8tt5pl2S71QnMFAKT+WS+qW3/zX6vk42NCQrTMp1206y9TnvTUP+KBs69v3IrFZgVCT447R9Oe8iC65z2z8VVrmOcsy/b/ask+Xs3xaE0zf0c5qYGZtUDyT1ufs/r/+gu+7f8Kk9eMuesc9p+P2MjbcCsZuCrkg4S6Nxgx9kIAQAzo9GBmVgAAABYhyAWAAAA1onz+XyB1zNjUHFxsdLS0jRQv1XjuMB7zqP+RPuyrFd1gqrMwv0ms1KBqWo6QiWz6L+5sYC5qUBVCUXuv0PNzQvM80vT3OV9c0l/01C3bVYkMK+5P82tbNDh2K1+9976mZvKkLjHPX7gpGKn3e4hN2Vh/XC3feTT7vK++f2Z35nX9yT5Vzrwqghh8nmkFnilmNS0KkZ9+sX3kz7Vu8oa+6rik5pFejgAwkBaQmCV8VpRUZFSU1M9+zETCwAAAOsQxAIAAMA6VCdAxISbQlCb6QdeVQiqc92Mf+922l5PwZvHE4zj5tP4P/R0n/43l/1zZ7tpBubmCJL/5gCbB7tVD/Zlu2kAFWn7jePuf/JNU930g9Ki5m67g1t1wKwi0Gi2e31JymzpjtFMhdjUwlz6ce+R/H3gfzObmzEk7DI3VHC/GzM9QvJPIfA/x+VVqcDrZ1ydFILa/D0C0PDUZi3ZaFTX6RLMxAIAAMA6BLEAAACwjnXpBI8++qjuvfdebdu2Td26ddPDDz+sE044IdLDQj2ozSXaUK7llQ6w4xz/wvvpzy5y2hUDega+2BGBl8bNpfRWy9yldzO1YN8At7pAQpH/ZVeNTzFelSkQs/JA7my3z5Y/u338l/rdKgK7Li6Wl5+LzXPcdtIG96hZYeBX8wKnPpjHzc0iTGbKhuSfhmG2/eoZGD9j8+dXmykApA0AqCtULjg8q2ZiZ86cqXHjxmnSpElatmyZunXrpiFDhuiHH/gfCQAAQENiVRD7wAMPaPTo0brssst0zDHHaNq0aWrWrJmeffbZSA8NAAAA9ciadIKysjItXbpUEyZMcI7Fx8dr0KBBWrRoUcBzSktLVVrqPsVcXOy9NApUZT6tbi5Bm+kDVYVSbD/DeJrevMePl/dz2mUt3HPNzQcyV7gVASRpS1rgVIHyie4S/a6v3E0JzDSFxp+6bfN+5gYHFcVuHzMtQZKaZu112j8b/cyND8yqBeamBklF7gYOZjqAuYmBV9UHyT9toOp7h0MKANDwsDQfm6yZid25c6fKy8vVunVrv+OtW7fWtm3bAp4zZcoUpaWlOX+ysrLqY6gAAACoY9bMxFbHhAkTNG7cOOd1cXExgSwAAA2MrfVYmUEOzpogNiMjQ40aNdL27dv9jm/fvl1t2rQJeE5SUpKSkpICvgc71OYGBzW5n5kOEOw/GrOfXzqC0aesS3unnWjcryTbrULQfoFbqcDc4KAox918QJKOmrLWaZtVE3Z/09ZpJxv9vdIGWn7l3vvfdz7htC/bcLI7jr908rv3unx3XBmfuAv8Jdlun6qbM1Qy0yK2/NlNgzjiRfc6ZpqA+Z1J/lUPZFQ0aL5mj9P22+DAYyOD6vx+eZ3TqOtRTrs6GycAAMJjTTpBYmKievXqpfnz5zvHKioqNH/+fPXr1y/ImQAAAIg11szEStK4ceM0cuRI9e7dWyeccIIefPBB7du3T5dddlmkhwYAAKIcy/Oxxaog9oILLtCOHTt02223adu2berevbvef//9Qx72AgAAQGyzKoiVpDFjxmjMmDGRHgbqSSTzYE1e5baqnmPmYZr9zFxZM5/T3RtLSvneLTuVYJThSvnezevefbxPJjMPtjTNzWtNyipxx2c+y/h1qtM082D3neP27/jBKPc6G9wRpuT53zvjE7ftPy7jcxjlwczSW2Zub7uH3HsU5bn9f/H4bJJ/ia7MOeucttfPz+tnIY/yXF67egW7R23lwdZ3HjgA2Mq6IBYAAKA6IlmlgFSG2mfNg10AAABAJWZiYY1oKbdV9d5maSWTV5qBuYuVqeqSeaUmP7qlsJK/9y+xtft4972EIvf8lDkpTnvnye79UvaY57pL8qOO/MJpz/77aca5bvmrlO/NfbL8r2umHZhpA62WuaXCzBQJc/eu9cPdc4982t2ly9zJq+pOaH7lz8zv3zjulRIQyrJ/JJfxSSEAgNAQxAIAANQBUgjqFukEAAAAsA4zsbBGXS2zhnLdYE+rlxnL3vELlwXsYy5hJxjL3+aT8uYT9+ZSurkLVdPZ/kvh267v77T3ZVcY7xiVCoylfnN5X8ZQZ39lphAETncwdxSTpLwZ7v3WD3fbiXvcfxuvy3fbOe3dezd5zO2fO9u9315j961f0t1zM1f5j8UrhcM8XuaRjhBKWgoVAgAg+hHEAgCABoHl/dhCOgEAAACsw0wsGqRwl4sPBHnq3e/J+SBpB5W80g+adO4bsL+5xF40uL/fe2mF5U671TJ3Wd7cUKEk261oYC7vp3zlphmY12n/ln8FhEqbhpb5vf6hZxOnXZHmvtdqmZsqkFTkbtSQ8Oz3TvvA8DbucaNqQcIu9/pJxvdUVZyxSUFJ/1ynnfLZeqedaPTx+vmZ6Qdmykewn3ddpBeQvgBUH7OrDRczsQAAALAOQSwAAACsQzoBGqRwUwhCPddredpcLjbTD8oG9HTaZhUCc7k8wW+jBP8l9k2nuUv6mStk9HNTCxL3mOkBgVMIvDZgMFMGMj7x/zevuRFCfJH7V0lRnlvFwNzAwauigFmJwT+1wG2XdWkf8FxJajp7sdP2edyjwjj/gJHC4fUXYH0v75NCgFjEMj/qGjOxAAAAsA5BLAAAAKxDOgFgMJe84zyebg/GTCGoMFIFDnhsgmCmCnjdr5HRp6rMFW5KgPlk/o5z8pz2r+a5j/yblQ5MZjUD04GTip12q4cS/d4ryW5StfvBa6W5bXNzBTNtwExfMCsmdL7M/f7Mn4VfBQj5pxf4VQ8wvkOTX6UC43h1UkaAhowUAUQTZmIBAABgHYJYAAAAWId0AsBgpgOEyutJdnMjA68+5hK7jLbfwr2xFB5fZXxJRsrCpgs7uceLfE57/flHOG0z/aAox61aYKYcmFKudu/93dWd/N4r7eBWJzA3SDCva1YYMFMZzHQC81wzBePHPHejhLR1/lUZzO/Wa4OJYBsWVPKqJmGqzmYHbF6AaMDSP2IdM7EAAACwDkEsAAAArEM6AVCLvJanvZa2fcZyu1d1ArN/1aVt8/y0QrfCwJYB7iYDLb9y+5vL+Enp7r9hzaV+c9MFswpAontYkpTyvXu/XReXOO3sSe49Ng8+QoGY9874926nvfOElk67yY8VTtvcdEGS2q8KnB5gCmVJP5QUguqkA5BCAAB1j5lYAAAAWIeZWAAAYlBOwZxDjvGwF2IJQSxQQ35L/8bxkDZOMDc1MPqb1zng0acqMw3gqM/c+62e0NFpp62T0S512uZmB+aSfkm2m5aQ+5p/BQOzskLJ16lOe2/nch2OmUJgXsdMITBTH9rMXux/ATMlw2PDCJNZ9SAhhBQO0gEAIPqRTgAAAADrMBMLAIAlSAcAXASxQC3yShXwKrbvd9xY2jaXyFWNpW2zqkDubHdTAnMp3VzGN6sZSO5GCe0Wusv7ZpqB5J+OEN+12GknLXC3ajCrEKQ/u8hpxxmfu9SjMoJfaoDH91eVV1WHxFWbAvYBANiLdAIAAABYh5lYAADqAakAQO2qVhC7adMmvfXWW9qwYYPKysr83nvggQdqZWBAQ2Muc1eYmwyYT+JXOcfraXwZqQLmRgGtlrldivKSnLa5IULmHLeEwaYLOzntX83zr05gbpCQPMe9VsIut/JAprGMv+n6/k47rdCtYPCLkXIg45pmO+Wz9QqF1yYFflUdQkgnqOlmB/V9XQBoiMIOYufPn6+hQ4eqY8eOWrVqlbp27arCwkL5fD717Nnz8BcAAAAAaijsIHbChAm64YYbNHnyZKWkpGjWrFlq1aqVLrroIp155pl1MUYAAKxXdfMB0guAmgk7iF25cqVefvnlgyc3bqyff/5ZzZs31x133KHf/va3uuqqq2p9kIDtvJ6a9zouYxneVP71ar/XjT2K/ptVCH41L3Bx/4xdgTcJKOmf67RbLfslYB9J+sWoVmBWKjCrHphpDWY6gnmtFOOa5mcw0xXMagtVxS9c5vleparfWyWvqhF1tdRPCgEA1J6wqxMkJyc7ebBt27bVunVu/tzOnTtrb2QAAACAh7BnYk888UR9+umnOvroo3X22Wdr/Pjx+s9//qM33nhDJ554Yl2MEQCAamHJHohdYQexDzzwgPbu3StJmjx5svbu3auZM2eqc+fOVCYAQuC1pFydpeZQNlfwekrfXPZPNJb3zQ0HzD5Vl/Qz/r07YD9zed/rfmaagZn64Jdm4FV5QVU2gzDu4ZU2YKJCAADEhrCD2I4dOzrt5ORkTZs2rVYHBABoeJgxBRCusHNiO3bsqF27dh1yfM+ePX4BLgAAAFBXwg5iCwsLVV5efsjx0tJSbd68uVYGBQAAAAQTcjrBW2+95bQ/+OADpaWlOa/Ly8s1f/585eTk1OrgABvUVo5lKNepWhLKzA0180nN/NUEjz7xRv5o2QB3o5LSlglO28yPrVpia8c5eU67yY8VTjvJuFaF0d8r99VvrEZ/r7Jfkv8OXqF8V2afUMqdof5VraFaU6QnALEv5CB22LBhkqS4uDiNHDnS772EhATl5OTo/vvvr9XBme666y7NmTNHK1asUGJiovbs2VNn9wIAAEB0CzmIrag4OKeSm5urL774QhkZGXU2qEDKysp0/vnnq1+/fnrmmWfq9d4AAACILmFXJ1i/3l3G++WXX9SkSZNaHZCXyZMnS5JmzJhRL/cDQlVby9ChXOeQPh5L42aZKzOD3exjtsuMPl5L9T8P7+t368w57kYnwXbUqrTT2OErc46bKpDosTuZyUxrqKrCSF8wP3coaQPhpiKg7rD8DyBcYT/YVVFRob/+9a/61a9+pebNm+u7776TJN16663MkAIAAKBehB3E3nnnnZoxY4buueceJSYmOse7du2qp59+ulYHV1OlpaUqLi72+wMAAAD7hZ1O8Pzzz+vJJ5/U6aefriuvvNI53q1bN61atSqsaxUUFGjq1KlB+6xcuVJdunQJd5iSpClTpjhpCECs81o+N5nVDMzdrcwKBmZqQKlRFaDqkr5fP4+KBns7t3Da5g5fZrWBpN37nbaZWuC3K1cV5mc1x+6VOkFKQPQLpzoBqQcApGoEsZs3b1anTp0OOV5RUaH9+/cHOMPb+PHjlZ+fH7RPTTZQmDBhgsaNG+e8Li4uVlZWVrWvBwAAgOgQdhB7zDHH6JNPPlF2drbf8ddff109evQI61qZmZnKzMwMdwghS0pKUlJSUp1dHwAAAJERdhB72223aeTIkdq8ebMqKir0xhtvaPXq1Xr++ef1zjvv1MUYJUkbNmzQ7t27tWHDBpWXl2vFihWSpE6dOql58+Z1dl8gmoWyZG4uy5vL+yk7jJQDY2OBUrOPUamgKnNjgoRDd6KW5F1VoOnsxU7brC7glxJh9A+WDmCmRdRWCgHpB7WL5X8AdSHsB7t++9vf6u2339a//vUvJScn67bbbtPKlSv19ttva/DgwXUxRkkHg+cePXpo0qRJ2rt3r3r06KEePXpoyZIldXZPAAAARKewZ2Il6eSTT9a8efNqeyxBzZgxgxqxABDlmHUFUF+qFcRKB3fQ+uGHH5ydvCp16NChxoMC4Aq2RB5KCkGckSrQ1Fh6l8d1m87+IWCfqszret3PPL7/iGSnndD1KKfd2KhIcMC4vt8GClU2UzA3NfD6fqhOAACxLewgds2aNbr88sv12Wef+R33+XyKi4tTeXm5x5kAAABA7Qg7iM3Pz1fjxo31zjvvqG3btoqLi6uLcQEALJRTMIeUAgD1IuwgdsWKFVq6dGm1NyAAEJ5gS+HmkrnXRgZe/b2W3s3r+IxzzcoGkn+FgUYeGxOYqQUy0gnM454pAB4pA5LkM9IRDnh8Vr/NEkgzAICYE3Z1gmOOOUY7d+6si7EAAAAAIQlpJra4uNhpT506VTfddJP+9re/6bjjjlNCQoJf39TU1NodIQDACqQRAKhPIQWxLVq08Mt99fl8Ov300/368GAXUP+8NgcIqb+xrG5WAkg0qwUY/ZvLXcKXpDgzBcE47vNILfihZxOn/atdbh+vTQ2CLvuHkAbglVJR0xQC0hEAIDqEFMR++OGHTruwsFBZWVlq1KiRX5+Kigpt2LChdkcHAAAABBBSEDtgwACnfdppp2nr1q1q1cr/QYtdu3Zp0KBBGjlyZO2OEAAAAKgi7OoElWkDVe3du1dNmjQJcAaA+uC1FO/FXPY3Nw84EKizqlQaqMLcyMC8VsWAnk47rfDwqUZeVRKq/kUVSiWGukIKAQBEh5CD2HHjxkmS4uLidOutt6pZs2bOe+Xl5Vq8eLG6d+9e6wMEAAAAqgo5iF2+fLmkgzOx//nPf5SYmOi8l5iYqG7duumGG26o/RECACKKqgMAolHIQWzlw12XXXaZ/vGPf1BKC4iAYGkCXsvc5tK7mRIQ7oYIVa/fyNhwIGHXPvceHhsWpBjHS/rnOu2mxjjMewdNE6jBkj7VBQAgNoSdEzt9+vS6GAcAIIKYbQVgm7B37AIAAAAiLeyZWACRE+ryt1+/MJfMzXPNlIGqy/vlHmkAoVw3abe7uYJXRYJQx+11b6/vKtwUgqrXJwUBAKIDQSwAQDkFc/xek14AINqRTgAAAADrMBMLNFChLMMH2+DA6yl/c4ODRKO/3+YFqzYFPNf8CynUSgxe6Q81qWDgdS8AQPQgiAUAOEgjAGAL0gkAAABgHWZigSgRySL8NX2S32/p3yNVoKKLW5Eg0UhTCCWtIZTqB1Lw9Idw7gcAiH4EsQDQQJE6AMBmpBMAAADAOgSxAAAAsA7pBECUiJacTK9duoLtXOWVz2tey6usltnHK6e1WjuV1aBPQ2FucEBqAQDbMBMLAAAA6zATCwANEDOvAGxHEAs0UF7L6mYKgd8OWFWW+r1KVfmlFnhdy2NXML9dvUIYayREshQaAMBFOgEAAACsQxALAA1QTsEcvwe7AMA2pBMA8OOVDhDyOV47bYWwS1ewaghe9zOFUjGhpkghAIDowEwsAAAArMNMLADEGCoPAGgICGKBGBRsKb5STZfeg1UuCHStUDZRqM44vPqx7A8AsY10AgAAAFiHIBYAYgyVBwA0BKQTADEo3KX0UPpXTVHwqlzglRJQ7tG/WvcmVQAAGjxmYgEAAGAdZmIBIEZVTSmgagGAWEIQCyAkoS7h+zLT3RdhpgrUtGJCXW1wAACIPlakExQWFmrUqFHKzc1V06ZNlZeXp0mTJqmsrCzSQwMAAEAEWDETu2rVKlVUVOiJJ55Qp06d9PXXX2v06NHat2+f7rvvvkgPDwAAAPXMiiD2zDPP1Jlnnum87tixo1avXq3HH3+cIBaoBV6bIwRbkg+3CkEoGzCE0icYUggAoOGwIp0gkKKiIrVs2TJon9LSUhUXF/v9AYCGitqxAGKJlUHs2rVr9fDDD+uKK64I2m/KlClKS0tz/mRlZdXTCAEAAFCXIppOUFBQoKlTpwbts3LlSnXp0sV5vXnzZp155pk6//zzNXr06KDnTpgwQePGjXNeFxcXE8gCAXhVBTBVZ6m/Udej3HuEmWYQbEykDQAAIhrEjh8/Xvn5+UH7dOzY0Wlv2bJFp556qvr3768nn3zysNdPSkpSUlJSTYcJADGD2rEAYkVEg9jMzExlZmaG1Hfz5s069dRT1atXL02fPl3x8VZmQgAAAKAWWFGdYPPmzRo4cKCys7N13333aceOHc57bdq0ieDIgNjjtVRf9XhIaQc7fjxsH6/7makIVa9TmxskAADsZEUQO2/ePK1du1Zr165V+/bt/d7z+XwRGhUA2I1UAgA2s2JNPj8/Xz6fL+AfAAAANDxWzMQCiG4VA3q6L1ZtcprhVj3w2jQhWHWCcK8LAIgNBLEA0ECQPgAglliRTgAAAACYmIkFUC3mEn1CZrrT9hntxh79zaX+kv65Tjvls8D9g6UDkCrgjZlXALGMmVgAAABYhyAWAAAA1iGdAEBIT/IHqxBgVhLwOserikDT2Yudts/c4KAaaQJUJPBnbjFLagGAWMNMLAAAAKxDEAsAAADrkE4AIKQUgtpcnve6VtyOHw87jmDnk0LgLadgDikFAGIKM7EAAACwDkEsAAAArEMQCwANhFmtAABsR04sAE+h5ph6lc8yhbKTV03HAQBoOJiJBQAAgHUIYgGggaA6AYBYQjoBgDpjphCYu3qZ6QcVA3o67QMLl9XPwAAA1mMmFgAAANYhiAUAAIB1SCcA4CnUHbvM9xp1PcppmztweVUwSNi1z2mXV2uUgdXVbmMAgOjATCwAAACsQxALAAAA65BOAMCPmQ5wwKgoECozhcArzcCvUkGYY6p6vhdSCAAgthHEAkCMoz4sgFhEOgEAAACsw0wsAD9eS/VVqwuEu1y//4hkp+31r2evlINQ0gcAAA0LQSwAxICvJw9RampqpIcBAPWGdAIAAABYh5lYAJ6qs2GA2c/v/IXLAvYv69LeaSeu2hSwT3WqEwAAYhtBLADEgK6TPlB8UjPP96lQACDWkE4AAAAA6zATC8BTTTcM8Eot2HFOntNOf3aR29/jOqQP1AyzsABiETOxAAAAsA5BLAAAAKxDOgGAemGmFqQ/W7M0BRweKQQAYh0zsQAAALAOQSwAAACsQzoBgFrltUFCKMfNjQ/iPTZHCHYtHEQqAYCGgJlYAAAAWIcgFgAAANaxJp1g6NChWrFihX744Qelp6dr0KBBmjp1qtq1axfpoQEw+DLT3RfGUr9XCoHJK4WgUdej/F4fYPODgEgjANCQWDMTe+qpp+rVV1/V6tWrNWvWLK1bt06/+93vIj0sAAAARIA1M7HXX3+9087OzlZBQYGGDRum/fv3KyEhIYIjAwAAQH2zJog17d69Wy+99JL69+8fNIAtLS1VaWmp87q4uLg+hgc0aOXGUr9X2oCZWuCXKuBRaaCc9IGQ5BTMCXicNAMAsciadAJJ+stf/qLk5GQdccQR2rBhg958882g/adMmaK0tDTnT1ZWVj2NFAAAAHUpokFsQUGB4uLigv5ZtWqV0//GG2/U8uXLNXfuXDVq1EiXXnqpfD6f5/UnTJigoqIi58/GjRvr42MBQEQV3n2O3x8AiEURTScYP3688vPzg/bp2LGj087IyFBGRoaOPPJIHX300crKytLnn3+ufv36BTw3KSlJSUlJtTlkAAF4bT4QykYEoaQKVK1OQHoBACCiQWxmZqYyMzOrdW5FRYUk+eW8AgAAoGGw4sGuxYsX64svvtBJJ52k9PR0rVu3Trfeeqvy8vI8Z2EBoKEhdQBAQ2LFg13NmjXTG2+8odNPP11HHXWURo0apeOPP14LFy4kXQAAAKABsmIm9rjjjtOCBQsiPQwAHrxyX71yZb2OeyEHFgBQlRUzsQCA4L6ePCTSQwCAekUQCwAAAOtYkU4AwE6+zHSnbf5l45VaYPYnhQAAEAxBLADEgK6TPlB8UjO/Y1QrABDLSCcAAACAdZiJBVBnzJQAv123PHb14i8kAECo+H8GAMSArycPUWpqaqSHAQD1hnQCAAAAWIeZWABRI5SNDxCY+WAXD3QBaAiYiQUAAIB1CGIBAABgHdIJANQLNi+oPzkFcw45RooBgFjDTCwAAACsQxALAAAA65BOACBq/Dy8r9NO+Wy906ZqweFRJxZAQ8NMLAAAAKxDEAsAAADrkE4AIGo0nb3YaR+I4DgAANGPmVgAAABYhyAWAAAA1iGdAEC9qxjQ02nHL1wWwZEAAGzFTCwAAACsQxALAAAA65BOAKDehZtCsGFyf7/XHSZ9VpvDAQBYiJlYAAAAWIcgFgAAANYhnQBA1CN9AABQFTOxAAAAsA5BLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxDEAsAAADrEMQCAADAOgSxAAAAsA5BLAAAAKxjXRBbWlqq7t27Ky4uTitWrIj0cAAAABAB1gWxN910k9q1axfpYQAAACCCrApi33vvPc2dO1f33XdfpIcCAACACGoc6QGEavv27Ro9erT+93//V82aNQvpnNLSUpWWljqvi4qKJEkHtF/y1ckwAaBe7VeZJKm4uDjCIwGA2lH595nPFzxYsyKI9fl8ys/P15VXXqnevXursLAwpPOmTJmiyZMnH3L8U71byyMEgMjKysqK9BAAoFaVlJQoLS3N8/043+HC3DpUUFCgqVOnBu2zcuVKzZ07V6+++qoWLlyoRo0aqbCwULm5uVq+fLm6d+/ueW7Vmdg9e/YoOztbGzZsCPqlwFVcXKysrCxt3LhRqampkR6ONfjewsd3Vj0//vijcnJyVFhYqPT09EgPBwBqzOfzqaSkRO3atVN8vHfma0SD2B07dmjXrl1B+3Ts2FG///3v9fbbbysuLs45Xl5erkaNGumiiy7Sc889F9L9iouLlZaWpqKiIv4nGSK+s+rhewsf31n18L0BaKgimk6QmZmpzMzMw/Z76KGHdOeddzqvt2zZoiFDhmjmzJnq27dvXQ4RAAAAUciKnNgOHTr4vW7evLkkKS8vT+3bt4/EkAAAABBBVpXYqqmkpCRNmjRJSUlJkR6KNfjOqofvLXx8Z9XD9wagoYpoTiwAAABQHQ1qJhYAAACxgSAWAAAA1iGIBQAAgHUIYgEAAGCdBh3EzpkzR3379lXTpk2Vnp6uYcOGRXpI1igtLVX37t0VFxenFStWRHo4UauwsFCjRo1Sbm6umjZtqry8PE2aNEllZWWRHlrUefTRR5WTk6MmTZqob9+++ve//x3pIUWtKVOmqE+fPkpJSVGrVq00bNgwrV69OtLDAoB61WCD2FmzZumSSy7RZZddpi+//FL/93//pwsvvDDSw7LGTTfdpHbt2kV6GFFv1apVqqio0BNPPKFvvvlGf//73zVt2jTdfPPNkR5aVJk5c6bGjRunSZMmadmyZerWrZuGDBmiH374IdJDi0oLFy7UNddco88//1zz5s3T/v37dcYZZ2jfvn2RHhoA1JsGWWLrwIEDysnJ0eTJkzVq1KhID8c67733nsaNG6dZs2bp2GOP1fLly9W9e/dID8sa9957rx5//HF99913kR5K1Ojbt6/69OmjRx55RJJUUVGhrKwsXXvttSooKIjw6KLfjh071KpVKy1cuFCnnHJKpIcDAPWiQc7ELlu2TJs3b1Z8fLx69Oihtm3b6qyzztLXX38d6aFFve3bt2v06NF64YUX1KxZs0gPx0pFRUVq2bJlpIcRNcrKyrR06VINGjTIORYfH69BgwZp0aJFERyZPYqKiiSJ3ysADUqDDGIrZ8Buv/12TZw4Ue+8847S09M1cOBA7d69O8Kji14+n0/5+fm68sor1bt370gPx0pr167Vww8/rCuuuCLSQ4kaO3fuVHl5uVq3bu13vHXr1tq2bVuERmWPiooKjR07Vr/+9a/VtWvXSA8HAOpNTAWxBQUFiouLC/qnMkdRkm655Radd9556tWrl6ZPn664uDi99tprEf4U9S/U7+3hhx9WSUmJJkyYEOkhR1yo35lp8+bNOvPMM3X++edr9OjRERo5Ys0111yjr7/+Wq+88kqkhwIA9apxpAdQm8aPH6/8/PygfTp27KitW7dKko455hjneFJSkjp27KgNGzbU5RCjUqjf24IFC7Ro0aJD9mjv3bu3LrroIj333HN1OMroEup3VmnLli069dRT1b9/fz355JN1PDq7ZGRkqFGjRtq+fbvf8e3bt6tNmzYRGpUdxowZo3feeUcff/yx2rdvH+nhAEC9iqkgNjMzU5mZmYft16tXLyUlJWn16tU66aSTJEn79+9XYWGhsrOz63qYUSfU7+2hhx7SnXfe6bzesmWLhgwZopkzZ6pv3751OcSoE+p3Jh2cgT311FOdGf/4+JhaAKmxxMRE9erVS/Pnz3fK3FVUVGj+/PkaM2ZMZAcXpXw+n6699lrNnj1bH330kXJzcyM9JACodzEVxIYqNTVVV155pSZNmqSsrCxlZ2fr3nvvlSSdf/75ER5d9OrQoYPf6+bNm0uS8vLymAXysHnzZg0cOFDZ2dm67777tGPHDuc9Zhld48aN08iRI9W7d2+dcMIJevDBB7Vv3z5ddtllkR5aVLrmmmv0z3/+U2+++aZSUlKc3OG0tDQ1bdo0wqMDgPrRIINY6WCZo8aNG+uSSy7Rzz//rL59+2rBggVKT0+P9NAQQ+bNm6e1a9dq7dq1hwT6DbC6nacLLrhAO3bs0G233aZt27ape/fuev/99w952AsHPf7445KkgQMH+h2fPn36YdNcACBWNMg6sQAAALAbyXkAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLEAAACwDkEsAAAArEMQCwAAAOsQxAIAAMA6BLFALdi6dasuvPBCHXnkkYqPj9fYsWMD9nvttdfUpUsXNWnSRMcdd5zefffd+h0oAAAxgiAWqAWlpaXKzMzUxIkT1a1bt4B9PvvsM/3hD3/QqFGjtHz5cg0bNkzDhg3T119/Xc+jBQDAfmw7C4Rgx44dOu644/TnP/9ZN998s6SDQenAgQP13nvv6fTTT3f6Dhw4UN27d9eDDz7od40LLrhA+/bt0zvvvOMcO/HEE9W9e3dNmzatXj4HAACxgplYIASZmZl69tlndfvtt2vJkiUqKSnRJZdcojFjxvgFsMEsWrRIgwYN8js2ZMgQLVq0qC6GDABATGsc6QEAtjj77LM1evRoXXTRRerdu7eSk5M1ZcqUkM/ftm2bWrdu7XesdevW2rZtW20PFQCAmMdMLBCG++67TwcOHNBrr72ml156SUlJSZEeEgAADRJBLBCGdevWacuWLaqoqFBhYWFY57Zp00bbt2/3O7Z9+3a1adOmFkcIAEDDQBALhKisrEwXX3yxLrjgAv31r3/VH//4R/3www8hn9+vXz/Nnz/f79i8efPUr1+/2h4qAAAxj5xYIES33HKLioqK9NBDD6l58+Z69913dfnllzvVBlasWCFJ2rt3r3bs2KEVK1YoMTFRxxxzjCTpuuuu04ABA3T//ffrnHPO0SuvvKIlS5boySefjNRHAgDAWpTYAkLw0UcfafDgwfrwww910kknSZIKCwvVrVs33X333brqqqsUFxd3yHnZ2dl+aQevvfaaJk6cqMLCQnXu3Fn33HOPzj777Pr6GAAAxAyCWAAAAFiHnFgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGAdglgAAABYhyAWAAAA1iGIBQAAgHUIYgEAAGCd/wewcy38XBhpsAAAAABJRU5ErkJggg==", + "text/plain": [ + "Figure(PyObject
)" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "fig, axis = plt.subplots(2, 2, figsize=(8,8))\n", + "fig.suptitle(\"NUTS - 21-D Neal's Funnel\", fontsize=16)\n", + "\n", + "fig.delaxes(axis[1,2])\n", + "fig.subplots_adjust(hspace=0)\n", + "fig.subplots_adjust(wspace=0)\n", + "\n", + "axis[1,1].hist(x10_nuts, bins=100, range=[-6,2])\n", + "axis[1,1].set_yticks([])\n", + "\n", + "axis[2,2].hist(theta_nuts, bins=100, orientation=\"horizontal\", range=[-4, 2])\n", + "axis[2,2].set_xticks([])\n", + "axis[2,2].set_yticks([])\n", + "\n", + "axis[2,1].hist2d(x10_nuts, theta_nuts, bins=100, range=[[-6,2],[-4, 2]])\n", + "axis[2,1].set_xlabel(\"x10\")\n", + "axis[2,1].set_ylabel(\"theta\");" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6a5d0a38", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.9.0", + "language": "julia", + "name": "julia-1.9" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.9.0" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/src/inference/Inference.jl b/src/inference/Inference.jl index 019a9bb55..a853629fa 100644 --- a/src/inference/Inference.jl +++ b/src/inference/Inference.jl @@ -149,7 +149,7 @@ end function AbstractMCMC.sample( rng::AbstractRNG, model::AbstractModel, - sampler::Sampler{<:InferenceAlgorithm}, + sampler::AbstractSampler, N::Integer; chain_type=MCMCChains.Chains, resume_from=nothing, @@ -210,7 +210,7 @@ end function AbstractMCMC.sample( rng::AbstractRNG, model::AbstractModel, - sampler::Sampler{<:InferenceAlgorithm}, + sampler::AbstractSampler, ensemble::AbstractMCMC.AbstractMCMCEnsemble, N::Integer, n_chains::Integer;