PDA_fit_ou_model.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T01:54:19.851981Z",
     "start_time": "2022-05-19T01:54:19.836361Z"
    }
   },
   "source": [
    "# Fitting diffusion model without analytical solution"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#  Model fitting:\n",
    "\n",
    "\n",
    "***\n",
    "## Models with likelihood function：\n",
    "### -MLE:  $ \\textit{max}\\ L(\\theta|Y) $\n",
    "### -Bayesian estimation:$ P(\\theta|Y)\\propto L(\\theta|Y)P(\\theta)$\n",
    "\n",
    "\n",
    "***\n",
    "## Models withoout likelihood function:\n",
    "### Replacing the likelihood function with summary statistics\n",
    "### Simulate data: X, Observed data: Y\n",
    "### -Point estimation: $ \\textit{min}\\ S(X,Y) $\n",
    "### -Approximate bayesian computation: $ P(\\theta|Y)\\propto S(X, Y)P(\\theta)$\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "  \n",
    "\n",
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Commonly used summary statistics in DDM/SSM:\n",
    "### 1. Chi-square\n",
    "<img src=\"chi.jpg\" style=\"zoom:50%\">\n",
    "\n",
    "\n",
    "\n",
    "### 2. KS-test\n",
    "<img src=\"ks.jpg\" style=\"zoom:50%\">\n",
    "\n",
    "\n",
    "\n",
    "### 3. Quantile likelihood\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "### 4. CAF/WLS\n",
    "<img src=\"caf.jpg\" style=\"zoom:50%\">\n",
    "\n",
    "\n",
    "### 5. PDA(Probability density approximation)\n",
    "<img src=\"kde_fun.jpg\" style=\"zoom:50%\">\n",
    "$$ \\large{\\mathscr{L}(\\theta|x)\\ = \\ \\Sigma \\ f(x|X)}$$\n",
    "\n",
    "----"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:47:38.003158Z",
     "start_time": "2022-05-19T22:47:29.804892Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import numba\n",
    "import pandas as pd\n",
    "from numba import jit,njit\n",
    "from numba import prange\n",
    "from KDEpy import FFTKDE\n",
    "from scipy.optimize import differential_evolution as de"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:47:44.784710Z",
     "start_time": "2022-05-19T22:47:38.175039Z"
    }
   },
   "outputs": [],
   "source": [
    "import hddm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "----\n",
    "\n",
    "# Ornstein–Uhlenbeck process model：\n",
    "$ dx = (v + \\lambda X)dt + \\sigma dW$\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "<img src=\"ou.png\" style=\"zoom:60%\">\n",
    "\n",
    "# Leaky competing accumulator model:\n",
    "## Binary choice task:\n",
    "$ dx_{1} = (I_1 - k X_1 - \\beta X_2)dt + \\sigma dW$\n",
    "\n",
    "\n",
    "$ dx_{2} = (I_2 - k X_2 - \\beta X_1)dt + \\sigma dW$\n",
    "\n",
    "$ x_i= \\textit{max}\\ (0,\\ x_i)$\n",
    "\n",
    "\n",
    "<img src=\"lca.jpg\" style=\"zoom:50%\">\n",
    "\n",
    "# OU Model $\\approx$  LCA Model:\n",
    "\n",
    "$ dx = [(I_1-I_2) + (\\beta - k) X]dt + \\sigma dW$\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:48:30.573034Z",
     "start_time": "2022-05-19T22:48:30.541782Z"
    }
   },
   "outputs": [],
   "source": [
    "@jit(fastmath=True)\n",
    "def ou_single(v,a,z,g,ndt,max_steps):\n",
    "    n_steps = 0\n",
    "    x = a * z\n",
    "    b_u = a\n",
    "    b_l = 0\n",
    "    dt = 0.001\n",
    "    while (x > b_l and x < b_u and n_steps < max_steps):\n",
    "        dx = (v + g*x) * dt + np.sqrt(dt) * np.random.normal()\n",
    "        x += dx\n",
    "        n_steps += 1.0\n",
    "        t = dt * n_steps\n",
    "    rt = n_steps * dt\n",
    "    return rt + ndt if x >= b_l else -rt - ndt\n",
    "\n",
    "@njit(fastmath=True,parallel=True)\n",
    "def ou_sim(v,a,z,g,ndt,max_steps,nsim):\n",
    "    data = np.empty(nsim)\n",
    "    for n in prange(nsim):\n",
    "        data[n] = ou_single(v,a,z,g,ndt,max_steps)\n",
    "    return data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Steps of PDA:\n",
    "<img src=\"pda.png\" style=\"zoom:50%\">\n",
    "\n",
    "\n",
    "---\n",
    "\n",
    "# Silver bandwidth:\n",
    "<img src='silver.jpg' style=\"zoom:50%\">\n",
    "\n",
    "\n",
    "<img src='kde.gif' style=\"zoom:100%\">\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:48:50.608091Z",
     "start_time": "2022-05-19T22:48:50.467471Z"
    }
   },
   "outputs": [],
   "source": [
    "def log_like(para):\n",
    "    global sim_np\n",
    "    v,a,z,g,ndt = para\n",
    "    real_data = sim_np.copy()\n",
    "    \n",
    "    # classify data into correct and incorrect response\n",
    "    cor_logrt = real_data[real_data['acc']==1]['log_rt']\n",
    "    incor_logrt = real_data[real_data['acc']==0]['log_rt']\n",
    "    N = real_data.shape[0]\n",
    "    N1 = cor_logrt.shape[0]\n",
    "    N2 = incor_logrt.shape[0]\n",
    "    \n",
    "    # non-decision time must be smaller than minRT \n",
    "    if ndt > np.min(np.absolute(real_data['rt'])):\n",
    "        sum_llf = -np.inf\n",
    "    else:\n",
    "        # construct pesudo log-likelihood \n",
    "        sim_N = 2000\n",
    "        sim_data = ou_sim(v,a,z,g,ndt,15000,sim_N)\n",
    "        df = pd.DataFrame({'rt':np.absolute(sim_data),'acc':(np.sign(sim_data)+1)/2})\n",
    "        df['log_rt'] = 0.0\n",
    "        df.loc[df.rt > 0.0, 'log_rt'] = np.log(df.loc[df.rt > 0.0, 'rt'])\n",
    "        acc_df = list(df.groupby('acc'))\n",
    "        ## in case there is only one response in the simulated data\n",
    "        if len(acc_df) == 1:\n",
    "            x,y = FFTKDE(kernel='epa').fit(acc_df[0][1]['log_rt'].to_numpy()).evaluate(2048)\n",
    "            if acc_df[0][1]['acc'].sample().to_list()[0] == 0:\n",
    "                sum_llf = np.sum(np.log(np.interp(incor_logrt,x,y)*(N2/N)))\n",
    "            elif acc_df[0][1]['acc'].sample().to_list()[0] == 1:\n",
    "                sum_llf = np.sum(np.log(np.interp(cor_logrt,x,y)*(N1/N)))\n",
    "        else:\n",
    "            sim_N1 = acc_df[1][1]['log_rt'].to_numpy().shape[0]\n",
    "            sim_N2 = acc_df[0][1]['log_rt'].to_numpy().shape[0]\n",
    "            x_ic,y_ic = FFTKDE(kernel='epa',bw='silverman').fit(acc_df[0][1]['log_rt'].to_numpy()).evaluate(2048)\n",
    "            x_c,y_c = FFTKDE(kernel='epa',bw='silverman').fit(acc_df[1][1]['log_rt'].to_numpy()).evaluate(2048)\n",
    "            cor_llf = np.log(np.interp(cor_logrt,x_c,y_c))\n",
    "            incor_llf = np.log(np.interp(incor_logrt,x_ic,y_ic))\n",
    "            ## p(rt,choice | theta) = p(rt_correct | theta) * p(rt_incorrect | theta) * p(choice \\ theta)\n",
    "            sum_llf  = np.sum(cor_llf) + np.sum(incor_logrt) + np.log(np.power((sim_N1/sim_N),N1) * np.power((sim_N2/sim_N),N2))\n",
    "    \n",
    "    return -sum_llf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:04.644059Z",
     "start_time": "2022-05-19T22:49:04.628436Z"
    }
   },
   "outputs": [],
   "source": [
    "## simulate data, trial number = 300\n",
    "v = 0.5\n",
    "a = 1.9\n",
    "z = 0.5\n",
    "g = -0.2\n",
    "ndt = 0.4\n",
    "max_steps = 15000\n",
    "nsim = 300\n",
    "\n",
    "sim_data = ou_sim(v,a,z,g,ndt,max_steps,nsim)\n",
    "sim_df = pd.DataFrame({'rt':np.absolute(sim_data),'log_rt':np.log(np.absolute(sim_data)),'acc':(np.sign(sim_data)+1)/2})\n",
    "sim_np = sim_df.to_records()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:06.198499Z",
     "start_time": "2022-05-19T22:49:06.167254Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "280.6078536754261"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "log_like([v,a,z,g,ndt])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:24.699307Z",
     "start_time": "2022-05-19T22:49:12.712943Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\KLIPS\\anaconda3\\lib\\site-packages\\KDEpy\\bw_selection.py:282: UserWarning: Silverman's rule failed. Too many idential values. Setting bw = 0.00026037783132902033\n",
      "  warnings.warn(\n",
      "C:\\Users\\KLIPS\\anaconda3\\lib\\site-packages\\KDEpy\\bw_selection.py:282: UserWarning: Silverman's rule failed. Too many idential values. Setting bw = 0.00024993262693512124\n",
      "  warnings.warn(\n"
     ]
    }
   ],
   "source": [
    "## fit data\n",
    "bound = [(0,2),(0,2),(0.3,0.7),(-1,1),(0.1,0.5)]\n",
    "result =de(func = log_like, bounds = bound)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-20T02:07:33.070697Z",
     "start_time": "2022-05-20T02:07:33.055080Z"
    }
   },
   "outputs": [],
   "source": [
    "result_df = pd.DataFrame(result['x']).T\n",
    "result_df.columns = ['v','a','z','g','ndt']\n",
    "result_df['log-likelihood'] = -result['fun']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-20T02:07:33.492572Z",
     "start_time": "2022-05-20T02:07:33.476953Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>v</th>\n",
       "      <th>a</th>\n",
       "      <th>z</th>\n",
       "      <th>g</th>\n",
       "      <th>ndt</th>\n",
       "      <th>log-likelihood</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.455289</td>\n",
       "      <td>1.70759</td>\n",
       "      <td>0.542364</td>\n",
       "      <td>-0.186065</td>\n",
       "      <td>0.493031</td>\n",
       "      <td>-272.662714</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          v        a         z         g       ndt  log-likelihood\n",
       "0  0.455289  1.70759  0.542364 -0.186065  0.493031     -272.662714"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "result_df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "# Deeping learning and ABC\n",
    "<img src='deep.png' style=\"zoom:50%\">\n",
    "&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp; Cramer et al., 2020\n",
    "\n",
    "## 1. Deep learning model approximates likelihood function:\n",
    "<img src='lan.png' style=\"zoom:100%\">\n",
    "&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp; Fengler et al., 2021\n",
    "\n",
    "## 2. Deep learning model approximates parameter posterior distribution:\n",
    "<img src='gan.png' style=\"zoom:100%\">\n",
    "\n",
    "&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp;&emsp; Ramesh et al., 2022"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:25.414936Z",
     "start_time": "2022-05-19T22:49:25.399312Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'doc': 'Model formulation is described in the documentation under LAN Extension.Meant for use with the extension.',\n",
       " 'params': ['v', 'a', 'z', 'g', 't'],\n",
       " 'params_trans': [0, 0, 1, 0, 0],\n",
       " 'params_std_upper': [1.5, 1.0, None, 1.0, 1.0],\n",
       " 'param_bounds': [[-2.0, 0.3, 0.2, -1.0, 0.001], [2.0, 2.0, 0.8, 1.0, 2]],\n",
       " 'boundary': <function hddm.simulators.boundary_functions.constant(t=0)>,\n",
       " 'params_default': [0.0, 1.0, 0.5, 0.0, 0.001],\n",
       " 'hddm_include': ['z', 'g'],\n",
       " 'choices': [-1, 1],\n",
       " 'slice_widths': {'v': 1.5,\n",
       "  'v_std': 0.1,\n",
       "  'a': 1,\n",
       "  'a_std': 0.1,\n",
       "  'z': 0.1,\n",
       "  'z_trans': 0.2,\n",
       "  't': 0.01,\n",
       "  't_std': 0.15,\n",
       "  'g': 0.1,\n",
       "  'g_trans': 0.2,\n",
       "  'g_std': 0.1}}"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Fitting simulated data with HDDM 0.9.5\n",
    "sim_df['response'] = sim_df['acc']\n",
    "sim_df['subj_idx'] = 0\n",
    "model = 'ornstein'\n",
    "n_samples = 2000\n",
    "hddm.model_config.model_config[model]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:25.774312Z",
     "start_time": "2022-05-19T22:49:25.758687Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['z', 'g']"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hddm.simulators.model_config[model]['hddm_include']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:49:26.196185Z",
     "start_time": "2022-05-19T22:49:26.118062Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'v': 1.5, 'v_std': 0.1, 'a': 1, 'a_std': 0.1, 'z': 0.1, 'z_trans': 0.2, 't': 0.01, 't_std': 0.15, 'g': 0.1, 'g_trans': 0.2, 'g_std': 0.1}\n",
      "Supplied model_config specifies params_std_upper for  z as  None.\n",
      "Changed to 10\n"
     ]
    }
   ],
   "source": [
    "hddmnn_model = hddm.HDDMnn(sim_df,\n",
    "                           informative = False,\n",
    "                           include = hddm.simulators.model_config[model]['hddm_include'],\n",
    "                           p_outlier = 0,\n",
    "                           w_outlier = 0.1,\n",
    "                           model = model,\n",
    "                          is_group_model = False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:50:05.219038Z",
     "start_time": "2022-05-19T22:49:26.539936Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " [-----------------100%-----------------] 2000 of 2000 complete in 38.6 sec"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<pymc.MCMC.MCMC at 0x18b64406d60>"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hddmnn_model.sample(n_samples, burn = 1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:50:05.594044Z",
     "start_time": "2022-05-19T22:50:05.562792Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>mean</th>\n",
       "      <th>std</th>\n",
       "      <th>2.5q</th>\n",
       "      <th>25q</th>\n",
       "      <th>50q</th>\n",
       "      <th>75q</th>\n",
       "      <th>97.5q</th>\n",
       "      <th>mc err</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>v</th>\n",
       "      <td>0.315275</td>\n",
       "      <td>0.0827841</td>\n",
       "      <td>0.160047</td>\n",
       "      <td>0.258748</td>\n",
       "      <td>0.315424</td>\n",
       "      <td>0.371004</td>\n",
       "      <td>0.477962</td>\n",
       "      <td>0.003733</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>a</th>\n",
       "      <td>0.960298</td>\n",
       "      <td>0.095541</td>\n",
       "      <td>0.823803</td>\n",
       "      <td>0.882141</td>\n",
       "      <td>0.948006</td>\n",
       "      <td>1.02417</td>\n",
       "      <td>1.17136</td>\n",
       "      <td>0.00850458</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>z</th>\n",
       "      <td>0.526762</td>\n",
       "      <td>0.0210816</td>\n",
       "      <td>0.481986</td>\n",
       "      <td>0.512725</td>\n",
       "      <td>0.527455</td>\n",
       "      <td>0.541161</td>\n",
       "      <td>0.565442</td>\n",
       "      <td>0.000857071</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>g</th>\n",
       "      <td>0.0170634</td>\n",
       "      <td>0.510281</td>\n",
       "      <td>-0.934204</td>\n",
       "      <td>-0.358197</td>\n",
       "      <td>-0.0360376</td>\n",
       "      <td>0.427953</td>\n",
       "      <td>0.914117</td>\n",
       "      <td>0.0448579</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t</th>\n",
       "      <td>0.424242</td>\n",
       "      <td>0.0248227</td>\n",
       "      <td>0.366175</td>\n",
       "      <td>0.409047</td>\n",
       "      <td>0.427204</td>\n",
       "      <td>0.44205</td>\n",
       "      <td>0.465009</td>\n",
       "      <td>0.001909</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "        mean        std      2.5q       25q        50q       75q     97.5q  \\\n",
       "v   0.315275  0.0827841  0.160047  0.258748   0.315424  0.371004  0.477962   \n",
       "a   0.960298   0.095541  0.823803  0.882141   0.948006   1.02417   1.17136   \n",
       "z   0.526762  0.0210816  0.481986  0.512725   0.527455  0.541161  0.565442   \n",
       "g  0.0170634   0.510281 -0.934204 -0.358197 -0.0360376  0.427953  0.914117   \n",
       "t   0.424242  0.0248227  0.366175  0.409047   0.427204   0.44205  0.465009   \n",
       "\n",
       "        mc err  \n",
       "v     0.003733  \n",
       "a   0.00850458  \n",
       "z  0.000857071  \n",
       "g    0.0448579  \n",
       "t     0.001909  "
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "hddmnn_model.gen_stats()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2022-05-19T22:50:05.953415Z",
     "start_time": "2022-05-19T22:50:05.937791Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-317.7683148094533"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(hddmnn_model.bic - 300*5)/2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.8.5"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}