diff --git a/notebooks/DIY-Figure1.ipynb b/notebooks/DIY-Figure1.ipynb index d6735c4..feeb197 100644 --- a/notebooks/DIY-Figure1.ipynb +++ b/notebooks/DIY-Figure1.ipynb @@ -7,7 +7,7 @@ "source": [ "## DIY Figure 1\n", "\n", - "In this notebook, we shown how to use the theoretical approximation to obtain accurate estimates of the angular error for the paleopole estimation for different choices of parameters. " + "In this notebook, we show how to use the theoretical approximation to obtain accurate estimates of the angular error for the paleopole estimation for different choices of parameters. A full simulation will be more expensive to run, but our simulations showed that theoretical approximations and numerical simulations differ at most $1\\%$ from each other when there are no outliers in the sample. " ] }, { @@ -20,7 +20,7 @@ "name": "stderr", "output_type": "stream", "text": [ - "/tmp/ipykernel_4124/2909491422.py:4: DeprecationWarning: Please use `gaussian_filter` from the `scipy.ndimage` namespace, the `scipy.ndimage.filters` namespace is deprecated.\n", + "/tmp/ipykernel_8973/2909491422.py:4: DeprecationWarning: Please use `gaussian_filter` from the `scipy.ndimage` namespace, the `scipy.ndimage.filters` namespace is deprecated.\n", " from scipy.ndimage.filters import gaussian_filter\n" ] } @@ -51,7 +51,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": 2, "id": "e4029aad-7ac2-48e9-989b-23e60b16a5ef", "metadata": {}, "outputs": [], @@ -65,7 +65,7 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": 3, "id": "779bead0-ee1c-4a77-8278-01872fd24d38", "metadata": {}, "outputs": [ @@ -145,17 +145,20 @@ "4 5 1 87.642494 8.676167" ] }, - "execution_count": 12, + "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ + "# Create a template data frames with all the combinations of number of sites and samples per site\n", "N_flat = np.tile(np.arange(1, N_max+1, 1), (n0_max,1)).ravel()\n", "n0_flat = np.tile(np.arange(1, n0_max+1, 1), (N_max,1)).T.ravel()\n", "\n", "df = pd.DataFrame({'N' : N_flat, 'n0' : n0_flat})\n", "\n", + "# Compute theoretical error in paleopole estimation\n", + "\n", "df[\"error_kappa_theoretical\"] = df.apply(lambda row: smp.kappa_theoretical(smp.Params(N=row.N, \n", " n0=row.n0, \n", " kappa_within_site=kappa_within_site,\n", @@ -181,7 +184,7 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": 4, "id": "701d20fc-c043-4cb2-96b3-da632087245c", "metadata": {}, "outputs": [], @@ -194,7 +197,7 @@ " idx = ((np.tile(a, (len(a0),1)).T - a0)**2).argmin(axis=1)\n", " return a0[idx].reshape(A.shape)\n", "\n", - "def plot_angle_error(df, df_in=None, save_plot=True):\n", + "def plot_angle_error(df):\n", " \n", " fig, axes = plt.subplots()\n", " fig.set_size_inches(14, 10)\n", @@ -248,19 +251,12 @@ " cmap='viridis', \n", " cbar_title=\"RMSE for pole estimation (degrees)\",\n", " color_max=14, \n", - " ticks=[1, 4, 7, 10, 15, 20, 30, 40])\n", - "\n", - " \n", - " if save_plot:\n", - " plt.savefig(\"Figure1.png\", dpi=300, format=\"png\", bbox_inches='tight')\n", - " plt.savefig(\"Figure1.pdf\", format=\"pdf\", bbox_inches='tight')\n", - " else: \n", - " plt.show()" + " ticks=[1, 4, 7, 10, 15, 20, 30, 40])" ] }, { "cell_type": "code", - "execution_count": 14, + "execution_count": 5, "id": "0d68f2bc-00d9-4457-b923-551ed76b1644", "metadata": {}, "outputs": [ @@ -276,12 +272,18 @@ } ], "source": [ - "df_pivot = df.pivot('N', 'n0', 'error_angle_theoretical')\n", + "df_pivot = df.pivot(index='N', columns='n0', values='error_angle_theoretical')\n", "\n", - "plot_angle_error(df_pivot, \n", - " df_pivot,\n", - " save_plot=True)" + "plot_angle_error(df_pivot)" ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "ccb48e2d-c60d-490f-a594-137f9fc382f5", + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { @@ -300,7 +302,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.11" + "version": "3.10.9" } }, "nbformat": 4,