From b4f24245ce9fd89572df941e5400488790ce2081 Mon Sep 17 00:00:00 2001
From: Berend Bouvy <b.n.bouvy@student.tudelft.nl>
Date: Fri, 20 Sep 2024 15:39:00 +0200
Subject: [PATCH] Added generalized function for seeing the effect of x_hat
---
content/GA_1_3/Analysis_Solution.ipynb | 89 +++++++++++++++++---------
content/GA_1_3/functions.py | 44 ++++++++++++-
2 files changed, 100 insertions(+), 33 deletions(-)
diff --git a/content/GA_1_3/Analysis_Solution.ipynb b/content/GA_1_3/Analysis_Solution.ipynb
index 21778555..fb13c882 100644
--- a/content/GA_1_3/Analysis_Solution.ipynb
+++ b/content/GA_1_3/Analysis_Solution.ipynb
@@ -31,7 +31,7 @@
},
{
"cell_type": "code",
- "execution_count": 101,
+ "execution_count": 1,
"id": "181ccfd5",
"metadata": {
"tags": []
@@ -40,10 +40,10 @@
{
"data": {
"text/plain": [
- "<Token var=<ContextVar name='format_options' default={'edgeitems': 3, 'threshold': 1000, 'floatmode': 'maxprec', 'precision': 8, 'suppress': False, 'linewidth': 75, 'nanstr': 'nan', 'infstr': 'inf', 'sign': '-', 'formatter': None, 'legacy': 9223372036854775807, 'override_repr': None} at 0x0000026E661E8810> at 0x0000026E1CAB5140>"
+ "<Token var=<ContextVar name='format_options' default={'edgeitems': 3, 'threshold': 1000, 'floatmode': 'maxprec', 'precision': 8, 'suppress': False, 'linewidth': 75, 'nanstr': 'nan', 'infstr': 'inf', 'sign': '-', 'formatter': None, 'legacy': 9223372036854775807, 'override_repr': None} at 0x0000026B5AF069D0> at 0x0000026B0DBBC4C0>"
]
},
- "execution_count": 101,
+ "execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
@@ -95,7 +95,7 @@
},
{
"cell_type": "code",
- "execution_count": 102,
+ "execution_count": 2,
"id": "8683b5ed",
"metadata": {},
"outputs": [
@@ -154,7 +154,7 @@
},
{
"cell_type": "code",
- "execution_count": 103,
+ "execution_count": 3,
"id": "41c56f43",
"metadata": {},
"outputs": [
@@ -202,7 +202,7 @@
},
{
"cell_type": "code",
- "execution_count": 104,
+ "execution_count": 4,
"id": "1f28eba3",
"metadata": {
"tags": []
@@ -251,7 +251,7 @@
},
{
"cell_type": "code",
- "execution_count": 105,
+ "execution_count": 5,
"id": "9f025cfc-4f89-492d-ac26-f5b6381d0c70",
"metadata": {
"tags": []
@@ -324,7 +324,7 @@
},
{
"cell_type": "code",
- "execution_count": 106,
+ "execution_count": 6,
"id": "71cf2133-37a6-4536-82a6-42a46b8a1c66",
"metadata": {
"tags": []
@@ -420,7 +420,7 @@
},
{
"cell_type": "code",
- "execution_count": 107,
+ "execution_count": 7,
"id": "f02ed4c4",
"metadata": {
"tags": []
@@ -442,7 +442,7 @@
},
{
"cell_type": "code",
- "execution_count": 108,
+ "execution_count": 8,
"id": "edf14892",
"metadata": {
"tags": []
@@ -480,7 +480,7 @@
},
{
"cell_type": "code",
- "execution_count": 109,
+ "execution_count": 9,
"id": "6cdfb46b-1324-4c2b-8148-5a6a102ede2e",
"metadata": {
"tags": []
@@ -572,7 +572,7 @@
},
{
"cell_type": "code",
- "execution_count": 110,
+ "execution_count": 10,
"id": "e868e488",
"metadata": {
"tags": []
@@ -665,7 +665,7 @@
},
{
"cell_type": "code",
- "execution_count": 111,
+ "execution_count": 11,
"id": "2c27b4f3",
"metadata": {},
"outputs": [],
@@ -727,7 +727,7 @@
},
{
"cell_type": "code",
- "execution_count": 112,
+ "execution_count": 12,
"id": "3a3eb1a1",
"metadata": {
"tags": []
@@ -770,7 +770,7 @@
},
{
"cell_type": "code",
- "execution_count": 113,
+ "execution_count": 13,
"id": "4bcd395d",
"metadata": {
"tags": []
@@ -876,7 +876,7 @@
},
{
"cell_type": "code",
- "execution_count": 114,
+ "execution_count": 14,
"id": "396ac3a5",
"metadata": {},
"outputs": [
@@ -961,7 +961,7 @@
},
{
"cell_type": "code",
- "execution_count": 115,
+ "execution_count": 15,
"id": "163acdb3",
"metadata": {
"colab": {
@@ -1011,7 +1011,7 @@
},
{
"cell_type": "code",
- "execution_count": 116,
+ "execution_count": 16,
"id": "5d583bd8",
"metadata": {
"tags": []
@@ -1109,7 +1109,7 @@
},
{
"cell_type": "code",
- "execution_count": 117,
+ "execution_count": 17,
"id": "13cffb4c",
"metadata": {},
"outputs": [],
@@ -1161,7 +1161,7 @@
},
{
"cell_type": "code",
- "execution_count": 118,
+ "execution_count": 18,
"id": "d85b1826",
"metadata": {
"tags": []
@@ -1236,7 +1236,7 @@
},
{
"cell_type": "code",
- "execution_count": 119,
+ "execution_count": 19,
"id": "4a592ac1",
"metadata": {
"tags": []
@@ -1344,7 +1344,7 @@
},
{
"cell_type": "code",
- "execution_count": 120,
+ "execution_count": 20,
"id": "835eefc8",
"metadata": {
"colab": {
@@ -1465,7 +1465,7 @@
},
{
"cell_type": "code",
- "execution_count": 121,
+ "execution_count": 21,
"id": "2711da12",
"metadata": {},
"outputs": [],
@@ -1501,7 +1501,7 @@
},
{
"cell_type": "code",
- "execution_count": 122,
+ "execution_count": 22,
"id": "d9a41ea5",
"metadata": {},
"outputs": [],
@@ -1524,7 +1524,7 @@
},
{
"cell_type": "code",
- "execution_count": 123,
+ "execution_count": 23,
"id": "b3bb808e",
"metadata": {},
"outputs": [
@@ -1611,7 +1611,7 @@
},
{
"cell_type": "code",
- "execution_count": 124,
+ "execution_count": 24,
"id": "ec7c8bef",
"metadata": {},
"outputs": [
@@ -1646,7 +1646,7 @@
},
{
"cell_type": "code",
- "execution_count": 125,
+ "execution_count": 25,
"id": "104d155d",
"metadata": {},
"outputs": [
@@ -1681,7 +1681,7 @@
},
{
"cell_type": "code",
- "execution_count": 126,
+ "execution_count": 26,
"id": "1dc93ce9",
"metadata": {},
"outputs": [
@@ -1741,7 +1741,7 @@
},
{
"cell_type": "code",
- "execution_count": 127,
+ "execution_count": 27,
"id": "113f3809",
"metadata": {},
"outputs": [
@@ -1783,14 +1783,14 @@
},
{
"cell_type": "code",
- "execution_count": 128,
+ "execution_count": 28,
"id": "aa877dd5",
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
- "model_id": "992db51bdea94fd2ad6417c9b7b27246",
+ "model_id": "af70788944c34efbac504e6738895a86",
"version_major": 2,
"version_minor": 0
},
@@ -1807,7 +1807,7 @@
"<function __main__.update_plot(x0, x1, x2)>"
]
},
- "execution_count": 128,
+ "execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
@@ -1846,6 +1846,31 @@
"interact(update_plot, x0=x0_slider, x1=x1_slider, x2=x2_slider)\n"
]
},
+ {
+ "cell_type": "code",
+ "execution_count": 29,
+ "id": "d8b31540",
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "application/vnd.jupyter.widget-view+json": {
+ "model_id": "67692528ebda4a6f9d74fd0125d2d79d",
+ "version_major": 2,
+ "version_minor": 0
+ },
+ "text/plain": [
+ "interactive(children=(FloatSlider(value=1.1806136803794327, description='xhat_0', max=15.122548861077993, min=…"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "xhat_slider_plot(model_gnss['A'], model_gnss['y'], model_gnss['times'], model_gnss['Sigma_Y'])"
+ ]
+ },
{
"cell_type": "markdown",
"id": "3203d779",
diff --git a/content/GA_1_3/functions.py b/content/GA_1_3/functions.py
index 0d030438..e33bc32b 100644
--- a/content/GA_1_3/functions.py
+++ b/content/GA_1_3/functions.py
@@ -2,6 +2,8 @@
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
+import ipywidgets as widgets
+from ipywidgets import interact
def plot_model(d, alt_model=None):
"""Time Series of observations with model and CI.
@@ -88,7 +90,7 @@ def plot_residual_histogram(d):
fig, ax = plt.subplots(figsize = (7,5))
- ax.hist(e_hat, bins = 40, density=True, edgecolor='black')
+ ax.hist(e_hat, density=True, edgecolor='black')
x = np.linspace(np.min(e_hat), np.max(e_hat), num=100);
ax.plot(x,
norm.pdf(x, loc=0.0, scale = np.std(e_hat)),
@@ -105,3 +107,43 @@ def plot_residual_histogram(d):
return fig, ax
+
+def xhat_slider_plot(A, y, t, Sigma_y=None):
+ """Interactive plot of the solution space for x_hat."""
+ n, k = A.shape
+
+ if Sigma_y is None:
+ Sigma_y = np.eye(n)
+ xhat = np.linalg.inv(A.T @ np.linalg.inv(Sigma_y) @ A) @ A.T @ np.linalg.inv(Sigma_y) @ y
+ Sigma_xhat = np.linalg.inv(A.T @ np.linalg.inv(Sigma_y) @ A)
+ std_xhat = np.sqrt(np.diag(Sigma_xhat))
+
+ sliders = {}
+ for i in range(k):
+ sliders[f'xhat_{i}'] = widgets.FloatSlider(value=xhat[i],
+ min=xhat[i] - 10*std_xhat[i],
+ max=xhat[i] + 10*std_xhat[i],
+ step=0.1*std_xhat[i],
+ description=f'xhat_{i}')
+
+ def update_plot(**kwargs):
+ xhat_values = np.array([kwargs[f'xhat_{i}'] for i in range(k)])
+ y_fit = A @ xhat_values
+ W = np.linalg.inv(Sigma_y)
+ ss_res = (y - y_fit).T @ W @ (y - y_fit)
+
+ plt.figure(figsize=(10, 5))
+ plt.plot(t, y, 'o', label='data')
+ plt.plot(t, y_fit, label='y_fit')
+ plt.title(f'Mean of squared residuals: {ss_res:.2f}')
+ plt.ylabel('y')
+ plt.xlabel('t')
+ plt.grid()
+ plt.legend()
+ plt.show()
+
+ interact(update_plot, **sliders)
+
+# Example usage
+# A, y, t should be defined before calling this function
+# xhat_slider_plot(A, y, t)
\ No newline at end of file
--
GitLab