examples/T2-hypothesis-testing.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np \n",
    "\n",
    "from pycircstat2 import Circular, load_data"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Hypothesis Testing\n",
    "\n",
    "### Table of Contents\n",
    "\n",
    "- [rayleigh_test](#the-rayleigh-test)\n",
    "- [V_test](#the-v-test)\n",
    "- [omnibus_test](#the-hodges-ajne-test)\n",
    "- [batschelet_test](#the-batschelet-test)\n",
    "- [chisquare_test](#chi-square-test)\n",
    "- [kuiper_test](#kuipers-test)\n",
    "- [watson_test](#watsons-one-sample-u2-test)\n",
    "- [raospacing_test](#raos-spacing-test)\n",
    "- [symmetry_test](#symmetry-test-around-the-median)\n",
    "- [one_sample_test](#one-sample-test)\n",
    "- [watson_williams_test](#watson-williams-test-for-two-multisample)\n",
    "- [watson_u2_test](#watsons-u2-test-for-two-multisample-with-or-without-ties)\n",
    "- [wheeler_watson_test](#wheeler-and-watson-two-sample-test)\n",
    "- [wallraff_test](#wallraffs-two-sample-test-for-angular-dispersion)\n",
    "\n",
    "### See also\n",
    "\n",
    "Chapter 27 of Zar (2010) contains many examples and step-by-step guide of how to compute most of circular hypothesis testing. We replicated all those examples and figures in notebook [`B2-Zar-2010`](https://nbviewer.org/github/circstat/pycircstat2/blob/main/examples/B2-Zar-2010.ipynb) with `pycircstats2`."
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Testing for Uniformity"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Rayleigh Test \n",
    "\n",
    "`rayleigh_test(alpha)` tests $H_{0}: \\rho=0$ vs. $H_{A}: \\rho \\neq 0$, where $\\rho$ is the population mean vector length. If the Rayleigh Test rejects $H_0$ ($p<0.05$), then the population is not a uniform circular distribution, or there is a mean population direction.\n",
    "\n",
    "**NOTE**: The Rayleigh Test assumes the data is <mark>unimodal</mark>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Rayleigh's Test of Uniformity\n",
      "-----------------------------\n",
      "H0: ρ = 0\n",
      "HA: ρ ≠ 0\n",
      "\n",
      "Test Statistics  (ρ | z-score): 0.82522 | 5.44787\n",
      "P-value: 0.00185 **\n",
      "Bootstrap P-value: 0.00680 **\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "RayleighTestResult(r=np.float64(0.8252177448200448), z=np.float64(5.4478746109270455), pval=np.float64(0.0018516375077209267), bootstrap_pval=0.0068)"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import rayleigh_test\n",
    "\n",
    "d1 = load_data('D1', source='zar')['θ'].values[:]\n",
    "c1 = Circular(data=d1)\n",
    "\n",
    "rayleigh_test(c1.alpha, B=9999, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The V Test\n",
    "\n",
    "`V_test(angle, alpha)` is a modified Rayleigh test that tests $H_{0}: \\rho=0$ vs. $H_{A}: \\rho \\neq 0$ and has a mean angle ($\\mu$)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Modified Rayleigh's Test of Uniformity\n",
      "--------------------------------------\n",
      "H0: ρ = 0\n",
      "HA: ρ ≠ 0 and μ = {angle:.5f} rad\n",
      "\n",
      "Test Statistics: 9.49761\n",
      "P-value: 0.00001 ***\n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import V_test\n",
    "\n",
    "d7 = load_data('D7', source='zar')['θ'].values[:]\n",
    "c7 = Circular(data=d7)\n",
    "\n",
    "V, u, pval = V_test(angle=np.deg2rad(90), alpha=c7.alpha, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Hodges-Ajne Test\n",
    "\n",
    "\n",
    "`omnibus_test(alpha)` tests $H_0$: uniform vs. $H_A$: not unifrom. Also called Ajne's A Test, or \"omnibus test\" because it works well for unimodal, bimodal, and multimoodal distributions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Hodges-Ajne (\"omnibus\") Test for Uniformity\n",
      "-------------------------------------------\n",
      "H0: uniform\n",
      "HA: not unifrom\n",
      "\n",
      "Test Statistics: 0.42752\n",
      "P-value: 0.00434 **\n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import omnibus_test\n",
    "\n",
    "d8 = load_data('D8', source='zar')['θ'].values[:]\n",
    "c8 = Circular(data=d8)\n",
    "\n",
    "A, pval = omnibus_test(c8.alpha, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### The Batschelet Test\n",
    "\n",
    "`batschelet_test(alpha)` is a modified Hodges-Ajne Test that tests $H_0$: uniform vs. $H_A$: not unifrom but concentrated around an angle θ."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Batschelet Test for Uniformity\n",
      "------------------------------\n",
      "H0: uniform\n",
      "HA: not unifrom but concentrated around θ = 0.78540 rad\n",
      "\n",
      "Test Statistics: 5\n",
      "P-value: 0.00661 **\n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import batschelet_test\n",
    "\n",
    "C, pval = batschelet_test(angle=np.deg2rad(45), alpha=c8.alpha, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Goodness-of-Fit Tests for Uniformity\n",
    "\n",
    "#### Chi-Square Test\n",
    "\n",
    "`chisquare_test(alpha)` tests the goodness of fit of a theoretical circular frequency distribution to an observed one. Here it is used to test whether the data in the population are distributed unifromly around the circle. This method is for <mark>grouped</mark> data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Chi-Square Test of Uniformity\n",
      "-----------------------------\n",
      "H0: uniform\n",
      "HA: not uniform\n",
      "\n",
      "Test Statistics (χ²): 66.54286\n",
      "P-value: 0.00000 ***\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "ChiSquareTestResult(chi2=np.float64(66.54285714285714), pval=np.float64(5.518107289173823e-10))"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import chisquare_test\n",
    "\n",
    "\n",
    "d2 = load_data(\"D2\", source=\"zar\")\n",
    "c2 = Circular(data=d2[\"θ\"].values[:], w=d2[\"w\"].values[:])\n",
    "\n",
    "chisquare_test(c2.w, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`kuiper_test(alpha)`, `watson_test(alpha)`, and `raospacing_test(alpha)` are Goodness-of-fit tests for ungrouped data. P-values for these tests are computed through simulation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pycircstat2.hypothesis import kuiper_test, watson_test, rao_spacing_test\n",
    "\n",
    "pigeon = np.array([20, 135, 145, 165, 170, 200, 300, 325, 335, 350, 350, 350, 355])\n",
    "c_pigeon = Circular(data=pigeon)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Kuiper's Test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Kuiper's Test of Circular Uniformity\n",
      "------------------------------------\n",
      "\n",
      "Test Statistic: 1.5047\n",
      "P-value = 0.1691 \n"
     ]
    }
   ],
   "source": [
    "V, pval = kuiper_test(c_pigeon.alpha, n_simulation=9999, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Watson's one-sample U2 Test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Watson's One-Sample U2 Test of Circular Uniformity\n",
      "--------------------------------------------------\n",
      "\n",
      "Test Statistic: 0.1361\n",
      "P-value = 0.1369 \n"
     ]
    }
   ],
   "source": [
    "U2, pval = watson_test(c_pigeon.alpha, n_simulation=9999, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Rao's Spacing Test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Rao's Spacing Test of Circular Uniformity\n",
      "-----------------------------------------\n",
      "\n",
      "Test Statistic: 2.8261\n",
      "P-value = 0.0766\n",
      "\n"
     ]
    }
   ],
   "source": [
    "U, pval = rao_spacing_test(c_pigeon.alpha, n_simulation=9999, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Testing for Symmetry"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Symmetry Test (around the median)\n",
    "\n",
    "`symmetry_test(alpha)` tests $H_0$: symmetrical around $\\theta$ vs. $H_A$: not symmetrical around $\\theta$, where $\\theta$ is the median of the population."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Symmetry Test\n",
      "------------------------------\n",
      "H0: symmetrical around median\n",
      "HA: not symmetrical around median\n",
      "\n",
      "Test Statistics: 14.50000\n",
      "P-value: 0.66406 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import symmetry_test\n",
    "\n",
    "d9 = load_data('D9', source='zar')['θ'].values[:]\n",
    "c9 = Circular(data=d9)\n",
    "\n",
    "statistics, pval = symmetry_test(alpha=c9.alpha, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Testing for the Mean Angle"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### One-Sample Test\n",
    "\n",
    "`one_sample_test(alpha)` tests $H_{0}: \\mu_a=\\mu_0$ vs. $H_{A}: \\mu_a \\neq \\mu_0$ ,where $\\mu_{a}$ is the population mean angle and $\\mu_{0}$ is a specified angle. This test is simply observing whether  $\\mu_{0}$ lies within the confidence interval for $\\mu_{a}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "One-Sample Test for the Mean Angle\n",
      "----------------------------------\n",
      "H0: μ = μ0\n",
      "HA: μ ≠ μ0 and μ0 = 1.57080 rad\n",
      "\n",
      "Failed to reject H0:\n",
      "μ0 = 1.57080 lies within the 95% CI of μ ([1.41993 1.86297])\n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import one_sample_test\n",
    "\n",
    "reject_or_not = one_sample_test(angle=np.deg2rad(90), alpha=c7.alpha, verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Two-Sample or Multisample Test"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Watson-Williams Test for Two-/Multisample\n",
    "\n",
    "`watson_williams_test(circs)` tests $H_0$: $\\mu_1 = \\mu_2 = ... = \\mu_n$ vs. $H_A$: $\\mu_1 \\neq \\mu_2 \\neq ... \\neq \\mu_n$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Watson-Williams Test for multiple samples\n",
      "---------------------------------------------\n",
      "H0: all samples are from populations with the same angle.\n",
      "HA: all samples are not from populations with the same angle.\n",
      "\n",
      "Test Statistics: 1.86524\n",
      "P-value: 0.18701 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import watson_williams_test\n",
    "\n",
    "data = load_data(\"D11\", source=\"zar\")\n",
    "s1 = Circular(data=data[data[\"sample\"] == 1][\"θ\"].values[:])\n",
    "s2 = Circular(data=data[data[\"sample\"] == 2][\"θ\"].values[:])\n",
    "s3 = Circular(data=data[data[\"sample\"] == 3][\"θ\"].values[:])\n",
    "\n",
    "F, pval = watson_williams_test(circs=[s1, s2, s3], verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Watson's U2 Test for Two-/multisample with or without Ties\n",
    "\n",
    "`watson_U2_test(circs)` tests $H_0$: $\\mu_1 = \\mu_2 = ... = \\mu_n$ vs. $H_A$: $\\mu_1 \\neq \\mu_2 \\neq ... \\neq \\mu_n$ for data with or without ties"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pycircstat2.hypothesis import watson_u2_test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Watson's U2 Test for two samples\n",
      "---------------------------------------------\n",
      "H0: The two samples are from populations with the same angle.\n",
      "HA: The two samples are not from populations with the same angle.\n",
      "\n",
      "Test Statistics: 0.14574\n",
      "P-value: 0.11261 \n"
     ]
    }
   ],
   "source": [
    "# without ties\n",
    "d = load_data(\"D12\", source=\"zar\")\n",
    "c0 = Circular(data=d[d[\"sample\"] == 1][\"θ\"].values[:])\n",
    "c1 = Circular(data=d[d[\"sample\"] == 2][\"θ\"].values[:])\n",
    "U2, pval = watson_u2_test(circs=[c0, c1], verbose=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Watson's U2 Test for two samples\n",
      "---------------------------------------------\n",
      "H0: The two samples are from populations with the same angle.\n",
      "HA: The two samples are not from populations with the same angle.\n",
      "\n",
      "Test Statistics: 0.06123\n",
      "P-value: 0.59716 \n"
     ]
    }
   ],
   "source": [
    "# with ties\n",
    "d = load_data(\"D13\", source=\"zar\")\n",
    "c0 = Circular(data=d[d[\"sample\"] == 1][\"θ\"].values[:], w=d[d[\"sample\"] == 1][\"w\"].values[:])\n",
    "c1 = Circular(data=d[d[\"sample\"] == 2][\"θ\"].values[:], w=d[d[\"sample\"] == 2][\"w\"].values[:])\n",
    "U2, pval = watson_u2_test(circs=[c0, c1], verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Wheeler and Watson Two-sample Test\n",
    "\n",
    "`wheeler_watson_test(circs)` tests $H_0$: $\\mu_1 = \\mu_2 = ... = \\mu_n$ vs. $H_A$: $\\mu_1 \\neq \\mu_2 \\neq ... \\neq \\mu_n$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Wheeler and Watson Two/Multi-Sample Test\n",
      "---------------------------------------------\n",
      "H0: All samples are from populations with the same angle.\n",
      "HA: All samples are not from populations with the same angle.\n",
      "\n",
      "Test Statistics: 3.67827\n",
      "P-value: 0.15895 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import wheeler_watson_test\n",
    "\n",
    "d = load_data(\"D12\", source=\"zar\")\n",
    "c0 = Circular(data=d[d[\"sample\"] == 1][\"θ\"].values[:])\n",
    "c1 = Circular(data=d[d[\"sample\"] == 2][\"θ\"].values[:])\n",
    "\n",
    "W, pval = wheeler_watson_test(circs=[c0, c1], verbose=True)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Wallraff's Two-sample Test for Angular Dispersion"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wallraff test of angular distances / dispersion\n",
      "-----------------------------------------------\n",
      "\n",
      "Test Statistics: 18.50000\n",
      "P-value: 0.77510 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import wallraff_test\n",
    "\n",
    "d = load_data(\"D14\", source=\"zar\")\n",
    "c0 = Circular(data=d[d[\"sex\"] == \"male\"][\"θ\"].values[:])\n",
    "c1 = Circular(data=d[d[\"sex\"] == \"female\"][\"θ\"].values[:])\n",
    "U, pval = wallraff_test(angle=np.deg2rad(135), circs=[c0, c1], verbose=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Wallraff test of angular distances / dispersion\n",
      "-----------------------------------------------\n",
      "\n",
      "Test Statistics: 13.00000\n",
      "P-value: 0.17524 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.utils import time2float\n",
    "\n",
    "d = load_data(\"D15\", source=\"zar\")\n",
    "c0 = Circular(data=time2float(d[d[\"sex\"] == \"male\"][\"time\"].values[:]))\n",
    "c1 = Circular(data=time2float(d[d[\"sex\"] == \"female\"][\"time\"].values[:]))\n",
    "U, pval = wallraff_test(angle=np.deg2rad(time2float(['7:55', '8:15'])), circs=[c0, c1], verbose=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# The Angular Randomization Test (ART) for homogeneity of two groups "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Angular Randomization Test (ART) for Homogeneity\n",
      "-------------------------------------------------\n",
      "H0: The two samples come from the same population.\n",
      "HA: The two samples do not come from the same population.\n",
      "\n",
      "Observed Test Statistic: 105.98051\n",
      "P-value: 0.57343 \n"
     ]
    }
   ],
   "source": [
    "from pycircstat2.hypothesis import angular_randomisation_test\n",
    "\n",
    "d = load_data(\"D13\", source=\"zar\")\n",
    "c0 = Circular(data=d[d[\"sample\"] == 1][\"θ\"].values[:])\n",
    "c1 = Circular(data=d[d[\"sample\"] == 2][\"θ\"].values[:])\n",
    "\n",
    "T, pval = angular_randomisation_test(circs=[c0, c1], n_simulation=1000, verbose=True)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Last updated: 2025-03-11 17:57:42CET\n",
      "\n",
      "Python implementation: CPython\n",
      "Python version       : 3.12.9\n",
      "IPython version      : 8.31.0\n",
      "\n",
      "numpy      : 2.2.3\n",
      "pycircstat2: 0.1.12\n",
      "\n",
      "Watermark: 2.5.0\n",
      "\n"
     ]
    }
   ],
   "source": [
    "%load_ext watermark\n",
    "%watermark --time --date --timezone --updated --python --iversions --watermark"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": ".venv",
   "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.12.9"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}