2_Sample_from_posterior.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f224cb62",
   "metadata": {},
   "outputs": [],
   "source": [
    "from fsbi.density_estimator import MakePosterior\n",
    "from fsbi.analyse import ComputeMetrics, default_x, condition\n",
    "from fsbi.utils import get_density_estim_data, resample_old_data\n",
    "import matplotlib.pyplot as plt\n",
    "import torch\n",
    "import numpy as np\n",
    "import h5py\n",
    "import argparse\n",
    "import yaml\n",
    "import pickle\n",
    "from time import time"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1031dc1f",
   "metadata": {},
   "source": [
    "#### 1/ Polynomial search space"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c83e05ed",
   "metadata": {},
   "outputs": [],
   "source": [
    "save_dir = \"data_synapsesbi/bg_IF_EEEIIEII_6pPol/\"\n",
    "runs_path = \"runs_synapsesbi/bg_IF_EEEIIEII_6pPol/\"\n",
    "metrics = [\"rate\"] #which metric(s) the posterior is on, see previous notebook.\n",
    "round_name = \"pi3_r5to10Hz\"\n",
    "\n",
    "\n",
    "with open(runs_path + \"posterior_\" + round_name + \".pkl\", \"rb\") as file:\n",
    "    mk_post = pickle.load(file)\n",
    "lower_lim=torch.tensor([0.01, 0.01, -2., -2., -2., -2.,\n",
    "            0.01, 0.01, -2., -2., -2., -2.,\n",
    "            0.01, 0.01, -2., -2., -2., -2.,\n",
    "            0.01, 0.01, -2., -2., -2., -2.])\n",
    "upper_lim=torch.tensor([.1, .1, 2., 2., 2., 2.,\n",
    "        .1, .1, 2., 2., 2., 2.,\n",
    "        .1, .1, 2., 2., 2., 2.,\n",
    "        .1, .1, 2., 2., 2., 2.])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8d9d1497",
   "metadata": {
    "heading_collapsed": true,
    "jp-MarkdownHeadingCollapsed": true
   },
   "source": [
    "#### 1 bis/ MLP search space"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6151dc40",
   "metadata": {
    "hidden": true
   },
   "outputs": [],
   "source": [
    "metrics = [\"rate\"] #which metrics is the posterior on\n",
    "round_name = \"pi3_r5to10Hz\"\n",
    "save_dir = \"../data_synapsesbi/bg_CVAIF_EEIE_T4wvceciMLP/\"\n",
    "runs_path = \"../runs_synapsesbi/bg_CVAIF_EEIE_T4wvceciMLP/\"\n",
    "with open(runs_path + \"posterior_\" + round_name + \".pkl\", \"rb\") as file:\n",
    "    mk_post = pickle.load(file)\n",
    "    \n",
    "lower_lim=torch.tensor([0., 0., -1., -1., -1., -1., #etaEE, etaIE, WpreEE, WpostEE, WpreIE, WpostIE\n",
    "            -1., -1., -1., -1., -1., -1.,\n",
    "            -1., -1., -1., -1., -1., -1.,\n",
    "            -1., -1., -1., -1.])\n",
    "upper_lim=torch.tensor([1., 1., 1., 1., 1., 1., #etaEE, etaIE, WpreEE, WpostEE, WpreIE, WpostIE\n",
    "            1., 1., 1., 1., 1., 1.,\n",
    "            1., 1., 1., 1., 1., 1.,\n",
    "            1., 1., 1., 1.])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18a4e2c5",
   "metadata": {},
   "source": [
    "#### 2/ Sample posterior"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "28e8f3ba",
   "metadata": {},
   "outputs": [],
   "source": [
    "#How many samples to generate\n",
    "post_num_samples = 10\n",
    "\n",
    "#What should their metric values be (here the only metric is the rate of excitation, and we want it to be between 5 and 10)\n",
    "default_x = torch.tensor([[np.random.random()*5+5] for i in range(post_num_samples)])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "831a4132",
   "metadata": {},
   "outputs": [],
   "source": [
    "#The bounds on the plasticity parameters\n",
    "bounds = {'low':lower_lim,\n",
    "          'high':upper_lim}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7e790eec",
   "metadata": {},
   "outputs": [],
   "source": [
    "#The sampling per say\n",
    "samples = mk_post.rsample(\n",
    "        default_x,\n",
    "        bounds)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "477ed1d5",
   "metadata": {},
   "outputs": [],
   "source": [
    "#Check that the samples are within the acceptable range\n",
    "in_prior = [ torch.all(torch.logical_and(lower_lim < samples[i,:], samples[i,:] < upper_lim)) for i in range(len(samples)) ]\n",
    "print(np.sum(in_prior)/len(in_prior)*100, \"% of samples are inside the original prior\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c51c4db5",
   "metadata": {},
   "outputs": [],
   "source": [
    "#Add the nuisance parameter: here only the input rate\n",
    "concat_tens = torch.cat( (samples, torch.rand( len(samples),1 )*10+5), dim=1) #between 5 and 15Hz"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c9402064",
   "metadata": {},
   "outputs": [],
   "source": [
    "## save the posterior samples\n",
    "np.savez(save_dir + \"post_samples_\" + round_name + \".npz\", thetas=concat_tens, default_x=default_x)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b63cb38",
   "metadata": {},
   "source": [
    "Great, you now have new rules to test!\n",
    "\n",
    "Go to Simulate_samples to simulate these rules with Auryn"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}