correlation_lower_upper_bounds.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": [],
      "authorship_tag": "ABX9TyMJEH8Q7zoyvB+bv+FBGl6M",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/microprediction/correlationbounds/blob/main/correlation_lower_upper_bounds.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Brute force correlation lower/upper bounds\n",
        "\n",
        "Laboriously creates a table that can be stored in a database.\n",
        "The usage is that any time we have an observation of correlation we can just join in the lower and upper bounds. Note that I've hardwired in the threshold parameter 0.025 here and that might be in the table too perhaps."
      ],
      "metadata": {
        "id": "NZDOKbfi1zSj"
      }
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "1yGxzVw0zskv",
        "outputId": "fde8f99c-4a47-424d-ffd1-7142e9df65f1"
      },
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "(0.16500000000000004, 0.4125000000000001)"
            ]
          },
          "metadata": {},
          "execution_count": 7
        }
      ],
      "source": [
        "import numpy as np\n",
        "import pandas as pd\n",
        "from scipy.stats import norm\n",
        "\n",
        "# Generate a correlated sample with true correlation rho\n",
        "def generate_correlated_sample(n, true_corr):\n",
        "    x = np.random.normal(size=n)\n",
        "    # Create y using the correlation structure\n",
        "    y = true_corr * x + np.sqrt(1 - true_corr**2) * np.random.normal(size=n)\n",
        "    return x, y\n",
        "\n",
        "# Function to compute the empirical correlation distribution for a given true correlation\n",
        "def compute_measured_correlations(true_corr, n, num_simulations=10000):\n",
        "    correlations = []\n",
        "    for _ in range(num_simulations):\n",
        "        x, y = generate_correlated_sample(n, true_corr)\n",
        "        measured_corr = np.corrcoef(x, y)[0, 1]\n",
        "        correlations.append(measured_corr)\n",
        "    return correlations\n",
        "\n",
        "# Function to compute lower and upper bounds using Bayes' rule\n",
        "def bayesian_correlation_bounds(measured_corr, n, true_corrs, num_simulations=10000):\n",
        "    posterior_probabilities = []\n",
        "\n",
        "    # For each true correlation, compute the likelihood of the measured correlation\n",
        "    for true_corr in true_corrs:\n",
        "        measured_corrs = compute_measured_correlations(true_corr, n, num_simulations)\n",
        "        likelihood = norm.pdf(measured_corr, loc=np.mean(measured_corrs), scale=np.std(measured_corrs))\n",
        "        posterior_probabilities.append(likelihood)\n",
        "\n",
        "    # Normalize the probabilities to get a posterior distribution\n",
        "    posterior_probabilities = np.array(posterior_probabilities)\n",
        "    posterior_probabilities /= posterior_probabilities.sum()\n",
        "\n",
        "    # Compute bounds (for simplicity, use percentiles)\n",
        "    lower_bound_idx = np.argmax(np.cumsum(posterior_probabilities) > 0.025)\n",
        "    upper_bound_idx = np.argmax(np.cumsum(posterior_probabilities) > 0.975)\n",
        "\n",
        "    lower_bound = true_corrs[lower_bound_idx]\n",
        "    upper_bound = true_corrs[upper_bound_idx]\n",
        "\n",
        "    return lower_bound, upper_bound\n",
        "\n",
        "\n",
        "bayesian_correlation_bounds(measured_corr=0.3, true_corrs=np.linspace(-0.99, 0.99, 101), n=250)"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "# Function to generate the correlation bounds table\n",
        "def generate_correlation_bounds_table(n, measured_correlations, true_corr_range=np.linspace(-0.99, 0.99, 101)):\n",
        "    data = []\n",
        "    for measured_corr in measured_correlations:\n",
        "        lower_bound, upper_bound = bayesian_correlation_bounds(measured_corr, n, true_corr_range)\n",
        "        data.append({\n",
        "            'count': n,\n",
        "            'measured_correlation': round(measured_corr, 3),\n",
        "            'lower_corr_bound': round(lower_bound, 3),\n",
        "            'upper_corr_bound': round(upper_bound, 3)\n",
        "        })\n",
        "\n",
        "    # Convert the data to a pandas DataFrame for better readability\n",
        "    df = pd.DataFrame(data)\n",
        "    return df\n",
        "\n",
        "bounds_table = generate_correlation_bounds_table(n=100, measured_correlations=np.linspace(-0.9,0.9,21))\n",
        "print(bounds_table)"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "f28A0LhU0HHP",
        "outputId": "444ac5a1-abc2-4462-d324-f851bace9c1b"
      },
      "execution_count": 10,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "    count  measured_correlation  lower_corr_bound  upper_corr_bound\n",
            "0     100                 -0.90            -0.931            -0.832\n",
            "1     100                 -0.81            -0.851            -0.713\n",
            "2     100                 -0.72            -0.792            -0.574\n",
            "3     100                 -0.63            -0.713            -0.475\n",
            "4     100                 -0.54            -0.653            -0.356\n",
            "5     100                 -0.45            -0.574            -0.257\n",
            "6     100                 -0.36            -0.495            -0.158\n",
            "7     100                 -0.27            -0.416            -0.079\n",
            "8     100                 -0.18            -0.356             0.020\n",
            "9     100                 -0.09            -0.277             0.099\n",
            "10    100                 -0.00            -0.198             0.178\n",
            "11    100                  0.09            -0.099             0.277\n",
            "12    100                  0.18            -0.020             0.356\n",
            "13    100                  0.27             0.079             0.436\n",
            "14    100                  0.36             0.158             0.495\n",
            "15    100                  0.45             0.257             0.574\n",
            "16    100                  0.54             0.356             0.653\n",
            "17    100                  0.63             0.475             0.713\n",
            "18    100                  0.72             0.574             0.792\n",
            "19    100                  0.81             0.713             0.851\n",
            "20    100                  0.90             0.832             0.931\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [],
      "metadata": {
        "id": "zVdYu4_Y1wi1"
      }
    },
    {
      "cell_type": "markdown",
      "source": [],
      "metadata": {
        "id": "oJsiyhDL1xKE"
      }
    }
  ]
}