{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "07aa6c24-cec2-491c-9db1-1e476f1cb584",
   "metadata": {},
   "source": [
    "# Ising model analyses\n",
    "\n",
    "_Copyright (C) 2023-, Joseph T. Lizier._\n",
    "_Distributed under GNU General Public License v3_\n",
    "\n",
    "This notebook template sets up reading in the data files for the Ising model simulations.\n",
    "\n",
    "First code block is just to start up the JVM. You may need to change your paths below (use the AutoAnalyser to show them for you):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "0eb38cde-2a28-4740-acf7-a76df9babf4d",
   "metadata": {},
   "outputs": [],
   "source": [
    "from jpype import *\n",
    "import numpy\n",
    "import sys\n",
    "# Our python data file readers are a bit of a hack, python users will do better on this:\n",
    "sys.path.append(\"../../../demos/python\")\n",
    "import readIntsFile\n",
    "\n",
    "if (not isJVMStarted()):\n",
    "    # Add JIDT jar library to the path\n",
    "    jarLocation = \"../../../infodynamics.jar\"\n",
    "    # Start the JVM (add the \"-Xmx\" option with say 1024M if you get crashes due to not enough memory space)\n",
    "    startJVM(getDefaultJVMPath(), \"-ea\", \"-Djava.class.path=\" + jarLocation, convertStrings=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8f3cc152-2695-4830-915b-ba1bb0306b69",
   "metadata": {},
   "source": [
    "Next load in a sample data file and check that we can plot the spins at one time step ok.\n",
    "\n",
    "Make sure that you have either downloaded and unzipped the [sample data](https://github.com/jlizier/jidt-cc-extras/raw/main/course/Module08-InfoTheoryAndSelfOrganisation/Ising/Ising-data.zip) or generated sample data. If generating yourself, then fix this range of temperatures that we have data files for:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "80112a52-e42a-44c2-abbd-5d612905fec5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 19216 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# We may want to loop over the data files later:\n",
    "tempsToExamine = numpy.array([2.00, 2.20, 2.26, 2.27, 2.28, 2.29, 2.34, 2.54]);\n",
    "\n",
    "# Check one data file - the way we reshape it into a grid may be useful for you to copy later:\n",
    "sampleTemp = tempsToExamine[3];\n",
    "filenameTemplate = './spins-%.2f.txt'\n",
    "sampleDataFilename = filenameTemplate % sampleTemp\n",
    "\n",
    "import os.path\n",
    "if (not os.path.isfile(sampleDataFilename)):\n",
    "    raise SystemExit(\"Ising model exported data files (e.g. %s) do not appear to exist in the current folder -- please fix this!\" % sampleDataFilename)\n",
    "\n",
    "dataRaw = readIntsFile.readIntsFile(sampleDataFilename)\n",
    "# As numpy array:\n",
    "data = numpy.array(dataRaw)\n",
    "import math\n",
    "modelSize = int(math.sqrt(data.shape[1]));\n",
    "timeSteps = data.shape[0];\n",
    "dataInLayout = numpy.zeros((timeSteps, modelSize, modelSize), dtype=int);\n",
    "for t in range(timeSteps):\n",
    "    # This reshapes in a transposed fashion to Matlab, but that doesn't matter for our purposes\n",
    "    dataInLayout[t,:,:] = numpy.reshape(data[t,:], (modelSize, modelSize));\n",
    "\n",
    "# Let's plot the last time step to check that it looks ok:\n",
    "import matplotlib.pyplot as plt\n",
    "plt.imshow( dataInLayout[timeSteps-1,:,:] )\n",
    "plt.title('Spins at t=%d for T=%.2f' % (timeSteps, sampleTemp));"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "26bd6e3c-3dfa-4e42-a384-158783eec2e2",
   "metadata": {},
   "source": [
    "# Analysis 1 -- MI versus spatial distance and temperature\n",
    "\n",
    "This analysis examines how MI between cells changes with spatial distance between the cells (explored via distance across the x-axis) -- at which temperatures does the system exhibit long-correlations?\n",
    "\n",
    "The script also examines how this changes with temperature, to explore the phase transition as well. Can you see the effect of the phase transition brought out here? (For finite-sized systems, the transition point is usually offset from the theoretical critical point to higher temperatures).\n",
    "\n",
    "The script may take several minutes to run:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "39fb6614-ed61-4ada-89ce-fcd5a4a7f23e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Done T=2.00\n",
      "Done T=2.20\n",
      "Done T=2.26\n",
      "Done T=2.27\n",
      "Done T=2.28\n",
      "Done T=2.29\n",
      "Done T=2.34\n",
      "Done T=2.54\n",
      "Computed MI for all distances at each temperature\n"
     ]
    }
   ],
   "source": [
    "# Compute MI vs x-separation for each temp:\n",
    "misVsXDist = numpy.zeros((tempsToExamine.size, math.ceil(modelSize/2)));\n",
    "for ti in range(tempsToExamine.size):\n",
    "    # 0. Load/prepare the data:\n",
    "    dataRaw = readIntsFile.readIntsFile(filenameTemplate % tempsToExamine[ti])\n",
    "    data = numpy.array(dataRaw)\n",
    "    timeSteps = data.shape[0];\n",
    "\n",
    "    # Let's reshape the data so it's easier to work with\n",
    "    dataInLayout = numpy.zeros((timeSteps, modelSize, modelSize), dtype=int);\n",
    "    for t in range(timeSteps):\n",
    "        # This reshapes in a transposed fashion to Matlab, but that doesn't matter for our purposes\n",
    "        dataInLayout[t,:,:] = numpy.reshape(data[t,:], (modelSize, modelSize));\n",
    "    \n",
    "    # 1. Construct the calculator:\n",
    "    calcClass = JPackage(\"infodynamics.measures.discrete\").MutualInformationCalculatorDiscrete\n",
    "    calc = calcClass(2, 2, 0)\n",
    "    # 2. No other properties to set for discrete calculators.\n",
    "    for xDiff in range(1,math.ceil(modelSize/2)+1):\n",
    "        # We'll calculate MI for cells separated by xDiff\n",
    "        # 3. Initialise the calculator for (re-)use:\n",
    "        calc.initialise();\n",
    "        for x in range(modelSize):\n",
    "            for y in range(modelSize):\n",
    "                # Column indices start from 0 in Python:\n",
    "                source = dataInLayout[:,x,y];\n",
    "                destination = dataInLayout[:,((x+xDiff) % modelSize),y];\n",
    "                # 4. Supply the sample data:\n",
    "                calc.addObservations(JArray(JInt, 1)(source.tolist()), JArray(JInt, 1)(destination.tolist()));\n",
    "        # 5. Compute the estimate:\n",
    "        # Note that we're placing the MI for xDiff in array index xDiff-1. Will correct this when plotting (below)\n",
    "        misVsXDist[ti, xDiff-1] = calc.computeAverageLocalOfObservations();\n",
    "    print('Done T=%.2f' % tempsToExamine[ti])\n",
    "print('Computed MI for all distances at each temperature')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "0feb0c6d-aed2-4e1c-a385-2e4398d0cfc7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x76edfd3410f0>]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "[encoded data, 62728 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 38384 characters]",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot MI vs x-separation for each temp:plt.figure();\n",
    "plt.figure();\n",
    "for ti in range(tempsToExamine.size):\n",
    "    # Add MI vs distance for this temperature:\n",
    "    plt.scatter(range(1,math.ceil(modelSize/2)+1), misVsXDist[ti,], marker='x', label='T=%.2f' % tempsToExamine[ti])\n",
    "plt.legend();\n",
    "plt.title('MI between cells as a function of x-separation');\n",
    "plt.ylabel('MI (bits)');\n",
    "plt.xlabel('x difference (cells)');\n",
    "\n",
    "# And plot MI vs temp for one sample distance:\n",
    "distance = 15;\n",
    "plt.figure();\n",
    "plt.title('MI between cells at x-separation=%d versus temperature' % distance);\n",
    "plt.ylabel('MI (bits)');\n",
    "plt.xlabel('temperature');\n",
    "plt.plot(tempsToExamine, misVsXDist[:,distance-1], marker='x')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5ade5f12-a804-48f0-8be1-5af0e33b67ce",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
