approach_BA.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T14:08:51.527354200Z",
     "start_time": "2023-08-16T14:08:19.499543600Z"
    },
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import pycalibration\n",
    "pycalibration.install()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T14:08:57.629151300Z",
     "start_time": "2023-08-16T14:08:51.528347100Z"
    },
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from pycalibration import ca\n",
    "import time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T15:55:31.890999300Z",
     "start_time": "2023-08-16T15:55:31.885542400Z"
    },
    "id": "DhSbZL5eoQUQ"
   },
   "outputs": [],
   "source": [
    "def SKCE(kernel, predictions, outcomes):\n",
    "\n",
    "    predictions = np.array(predictions)\n",
    "    outcomes = np.array(outcomes)\n",
    "    n = len(outcomes)\n",
    "\n",
    "    # Initialize SKCE\n",
    "    skce_matrix = np.zeros((n,n))\n",
    "\n",
    "    # For each pair of data points\n",
    "    for i in range(n):\n",
    "        for j in range(i+1,n):\n",
    "            # Compute the difference between unit vectors outcomes and predictions\n",
    "            diff_Yi_g_Xi = np.array([outcomes[i] - predictions[i]])\n",
    "            diff_Yj_g_Xj = np.array([outcomes[j] - predictions[j]])\n",
    "\n",
    "            # Compute the kernel matrix for the pair\n",
    "            k_g_Xi_g_Xj = kernel(predictions, predictions)\n",
    "\n",
    "            matrix_kernel = k_g_Xi_g_Xj[i,j] * np.eye(len(predictions[i]))\n",
    "\n",
    "            skce_matrix[i,j] = ((diff_Yi_g_Xi @ matrix_kernel) @ diff_Yj_g_Xj.T)[0,0]\n",
    "\n",
    "    skce = np.sum(skce_matrix)\n",
    "\n",
    "    \n",
    "    # Divide the SKCE by n over 2\n",
    "    prefactor = (n* (n-1)) // 2\n",
    "    skce /= prefactor\n",
    "\n",
    "    return skce\n",
    "\n",
    "\n",
    "def transform_outcomes(predictions, outcomes):\n",
    "    #creates vector for outcomes to match predictions shape\n",
    "    outcomes_minus_one = [x - 1 for x in outcomes]\n",
    "    transformed_outcomes = np.zeros_like(predictions)\n",
    "    transformed_outcomes[np.arange(len(predictions)), outcomes_minus_one] = 1\n",
    "    return transformed_outcomes\n",
    "\n",
    "# Use the rbf_kernel function from sklearn.metrics.pairwise\n",
    "def rbf_pw(X, Y=None, length_scale=1.0):\n",
    "    gamma = 1.0 / (2 * length_scale ** 2)\n",
    "    return RBF_PW(X, Y, gamma=gamma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def SKCE_torch(kernel, predictions, outcomes):\n",
    "       predictions = torch.tensor(predictions, dtype=torch.float64)\n",
    "       outcomes = torch.tensor(outcomes, dtype=torch.float64)\n",
    "       n = len(outcomes)\n",
    "       diff_Y_g_X = outcomes - predictions\n",
    "\n",
    "       # kernel matrix\n",
    "       k_g_Xi_g_Xj = kernel(predictions, predictions)\n",
    "\n",
    "       # Generate all pairwise combinations of indices\n",
    "       i_indices, j_indices = torch.triu_indices(n, n, 1)\n",
    "\n",
    "       # Create an identity matrix for matrix_kernel\n",
    "       identity_matrix = torch.eye(predictions.shape[1], dtype=torch.float32)\n",
    "\n",
    "       # Select the pairwise differences\n",
    "       diff_Yi_g_Xi = diff_Y_g_X[i_indices]\n",
    "       diff_Yj_g_Xj = diff_Y_g_X[j_indices]\n",
    "\n",
    "       # Compute the matrix kernel for all pairs\n",
    "       # k_g_Xi_g_Xj[i_indices, j_indices] Selects Upper Triangle Pairs\n",
    "       matrix_kernels = k_g_Xi_g_Xj[i_indices,j_indices].unsqueeze(-1).unsqueeze(-1) * identity_matrix\n",
    "       print(matrix_kernels.shape)\n",
    "       # Compute the SKCE matrix by first Yi_Xi @ Upper Triangle and then @ Yj_Xj\n",
    "       skce_values = (diff_Yi_g_Xi.unsqueeze(1) @ matrix_kernels @ diff_Yj_g_Xj.unsqueeze(-1)).reshape(-1)\n",
    "\n",
    "       # Sum the SKCE values\n",
    "       skce = torch.sum(skce_values)\n",
    "\n",
    "       # Normalize the SKCE value\n",
    "       prefactor = (n * (n - 1)) // 2\n",
    "       skce /= prefactor\n",
    "\n",
    "       return skce"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T15:55:32.941583200Z",
     "start_time": "2023-08-16T15:55:32.914890400Z"
    },
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from sklearn.gaussian_process.kernels import RBF as RBF_GP\n",
    "from sklearn.metrics.pairwise import rbf_kernel as RBF_PW\n",
    "\n",
    "# Instantiate the RBF kernel from sklearn.gaussian_process.kernels\n",
    "rbf_gp = RBF_GP(length_scale=1)\n",
    "\n",
    "# Use the rbf_kernel function from sklearn.metrics.pairwise\n",
    "def rbf_pw(X, Y=None, length_scale=1.0):\n",
    "    gamma = 1.0 / (2 * length_scale ** 2)\n",
    "    return RBF_PW(X, Y, gamma=gamma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T15:55:34.007036600Z",
     "start_time": "2023-08-16T15:55:33.910304200Z"
    },
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "vYPFfze2mMkW",
    "outputId": "f9b22c07-5700-48a0-d528-08f1b67f9124"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.33333333333333337"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rng = np.random.default_rng(1234)\n",
    "import random\n",
    "rng = np.random.default_rng(1234)\n",
    "predictions = [rng.dirichlet((3, 2, 5)) for _ in range(40)]\n",
    "outcomes = rng.integers(low=1, high=4, size=40)\n",
    "predictions = np.array([[0.1,0.9],[0.9,0.1]])\n",
    "outcomes = np.array([2,1])\n",
    "predictions = np.array([[0.5,0.5],[0.5,0.5]])\n",
    "outcomes = np.array([2,2])\n",
    "\n",
    "\n",
    "\n",
    "#2 Datenpukte 3 Klassen\n",
    "predictions = np.array([[1/3,1/3,1/3],[1/3,1/3,1/3]])\n",
    "outcomes = np.array([3,1])\n",
    "\n",
    "predictions = np.array([[0.5,0.5],[0.5,0.5],[0.5,0.5]])\n",
    "outcomes = np.array([2,1,1])\n",
    "#outcomes = predictions.copy()\n",
    "predictions = np.array([[0.25,0.75],[0.5,0.5],[0.75,0.25]])\n",
    "outcomes = np.array([2,2,1])\n",
    "predictions = np.array([[1/3,1/3,1/3],[1/3,1/3,1/3]])\n",
    "outcomes = np.array([3,1])\n",
    "\n",
    "\n",
    "#print(transform_outcomes(predictions, outcomes))\n",
    "\n",
    "# Now you can use these in your calculate_SKCE function\n",
    "SKCE(rbf_pw, predictions, transform_outcomes(predictions, outcomes))\n",
    "#print(SKCE(rbf_pw, predictions, outcomes))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T15:55:36.323408100Z",
     "start_time": "2023-08-16T15:55:36.297604500Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1., 1.],\n",
       "       [1., 1.]])"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rbf_pw(predictions, predictions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-16T15:55:39.291692500Z",
     "start_time": "2023-08-16T15:55:39.270026200Z"
    },
    "id": "iUIjnGPqo8bA"
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.3333333333333333"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#skce_julia = ca.SKCE(ca.RBFKernel(), unbiased=True)\n",
    "skce_julia = ca.SKCE(ca.tensor(ca.RBFKernel(), ca.WhiteKernel()), unbiased = True)\n",
    "skce_julia(ca.RowVecs(predictions), outcomes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-02T15:53:04.126108100Z",
     "start_time": "2023-08-02T15:53:04.111226400Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0 0 0 0 1 1 1 2 2 3]\n"
     ]
    }
   ],
   "source": [
    "#skce_julia = ca.SKCE(ca.RBFKernel(), unbiased=True)\n",
    "#skce_julia(predictions, outcomes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-02T15:42:40.329970600Z",
     "start_time": "2023-08-02T15:42:40.151410700Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "0\n",
      "1\n",
      "0\n",
      "2\n",
      "1\n",
      "2\n"
     ]
    },
    {
     "ename": "TypeError",
     "evalue": "unsupported operand type(s) for +: 'NoneType' and 'NoneType'",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[1;32mIn[24], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[43mSKCE_vector\u001b[49m\u001b[43m(\u001b[49m\u001b[43mrbf_pw\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mpredictions\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mtransform_outcomes\u001b[49m\u001b[43m(\u001b[49m\u001b[43mpredictions\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43moutcomes\u001b[49m\u001b[43m)\u001b[49m\u001b[43m)\u001b[49m\n",
      "Cell \u001b[1;32mIn[21], line 20\u001b[0m, in \u001b[0;36mSKCE_vector\u001b[1;34m(kernel, predictions, outcomes)\u001b[0m\n\u001b[0;32m     18\u001b[0m vhelper \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39mvectorize(SKCE_vectorized_helper)\n\u001b[0;32m     19\u001b[0m results \u001b[38;5;241m=\u001b[39m vhelper(listed_indices_triu)\n\u001b[1;32m---> 20\u001b[0m skce \u001b[38;5;241m=\u001b[39m \u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43msum\u001b[49m\u001b[43m(\u001b[49m\u001b[43mresults\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m     23\u001b[0m \u001b[38;5;66;03m# Divide the SKCE by n over 2\u001b[39;00m\n\u001b[0;32m     24\u001b[0m prefactor \u001b[38;5;241m=\u001b[39m (n\u001b[38;5;241m*\u001b[39m (n\u001b[38;5;241m-\u001b[39m\u001b[38;5;241m1\u001b[39m)) \u001b[38;5;241m/\u001b[39m\u001b[38;5;241m/\u001b[39m \u001b[38;5;241m2\u001b[39m\n",
      "File \u001b[1;32m<__array_function__ internals>:200\u001b[0m, in \u001b[0;36msum\u001b[1;34m(*args, **kwargs)\u001b[0m\n",
      "File \u001b[1;32m~\\miniconda3\\lib\\site-packages\\numpy\\core\\fromnumeric.py:2324\u001b[0m, in \u001b[0;36msum\u001b[1;34m(a, axis, dtype, out, keepdims, initial, where)\u001b[0m\n\u001b[0;32m   2321\u001b[0m         \u001b[38;5;28;01mreturn\u001b[39;00m out\n\u001b[0;32m   2322\u001b[0m     \u001b[38;5;28;01mreturn\u001b[39;00m res\n\u001b[1;32m-> 2324\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m \u001b[43m_wrapreduction\u001b[49m\u001b[43m(\u001b[49m\u001b[43ma\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43madd\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[38;5;124;43msum\u001b[39;49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43maxis\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mdtype\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mout\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mkeepdims\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mkeepdims\u001b[49m\u001b[43m,\u001b[49m\n\u001b[0;32m   2325\u001b[0m \u001b[43m                      \u001b[49m\u001b[43minitial\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43minitial\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mwhere\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mwhere\u001b[49m\u001b[43m)\u001b[49m\n",
      "File \u001b[1;32m~\\miniconda3\\lib\\site-packages\\numpy\\core\\fromnumeric.py:86\u001b[0m, in \u001b[0;36m_wrapreduction\u001b[1;34m(obj, ufunc, method, axis, dtype, out, **kwargs)\u001b[0m\n\u001b[0;32m     83\u001b[0m         \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[0;32m     84\u001b[0m             \u001b[38;5;28;01mreturn\u001b[39;00m reduction(axis\u001b[38;5;241m=\u001b[39maxis, out\u001b[38;5;241m=\u001b[39mout, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mpasskwargs)\n\u001b[1;32m---> 86\u001b[0m \u001b[38;5;28;01mreturn\u001b[39;00m ufunc\u001b[38;5;241m.\u001b[39mreduce(obj, axis, dtype, out, \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m*\u001b[39mpasskwargs)\n",
      "\u001b[1;31mTypeError\u001b[0m: unsupported operand type(s) for +: 'NoneType' and 'NoneType'"
     ]
    }
   ],
   "source": [
    "# Measure the time taken by function1\n",
    "start_time = time.time()\n",
    "SKCE(rbf_gp, predictions, outcomes)\n",
    "end_time = time.time()\n",
    "print(f\"SKCE Python took {end_time - start_time} seconds to execute.\")\n",
    "\n",
    "# Measure the time taken by function2\n",
    "start_time = time.time()\n",
    "skce_julia(ca.RowVecs(predictions), outcomes)\n",
    "end_time = time.time()\n",
    "print(f\"SKCE Julia took {end_time - start_time} seconds to execute.\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 93,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-02T15:34:03.514720100Z",
     "start_time": "2023-08-02T15:34:03.462145100Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.03245003262797521"
      ]
     },
     "execution_count": 93,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from julia import Main\n",
    "import numpy as np\n",
    "\n",
    "def test_rbf_kernel(predictions):\n",
    "    Main.eval(\"\"\"\n",
    "    using KernelFunctions\n",
    "    rbf = SqExponentialKernel()\n",
    "    \"\"\")\n",
    "\n",
    "    kernel_matrix = np.zeros((len(predictions), len(predictions)))\n",
    "\n",
    "    for i in range(len(predictions)):\n",
    "        for j in range(i+1,len(predictions)):\n",
    "            kernel_matrix[i][j] = Main.eval(f\"rbf({list(predictions[i])}, {list(predictions[j])})\")\n",
    "\n",
    "    return kernel_matrix\n",
    "\n",
    "# Test the function with two predictions\n",
    "#predictions = np.array([[0.1,0.9],[0.9,0.1],[0.5,0.5]])\n",
    "print(test_rbf_kernel(predictions))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 94,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-08-02T15:34:04.120705400Z",
     "start_time": "2023-08-02T15:34:04.074023400Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.03245003262797521"
      ]
     },
     "execution_count": 94,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rbf_pw(predictions, predictions)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-20T13:24:50.151910500Z",
     "start_time": "2023-07-20T13:24:49.986245900Z"
    },
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.9394130628134758"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Main.eval(f\"rbf({list(predictions[0])}, {list(predictions[1])})\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-20T13:24:50.338902800Z",
     "start_time": "2023-07-20T13:24:50.180678400Z"
    },
    "collapsed": false
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "colab": {
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}