Skip to content

Commit

Permalink
kalman
Browse files Browse the repository at this point in the history
  • Loading branch information
littlehome authored and littlehome committed Oct 30, 2023
1 parent 8047536 commit 78a98a9
Show file tree
Hide file tree
Showing 3 changed files with 351 additions and 0 deletions.
168 changes: 168 additions & 0 deletions content/finance/2023-10-29-physics.md
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,168 @@
+++
title = "Physics"
author = ["littlehome"]
date = 2023-10-29
lastmod = 2023-10-30T15:43:45+09:00
draft = false
+++

## Fundamentals of physics {#fundamentals-of-physics}

<https://www.youtube.com/playlist?list=PLFE3074A4CB751B2B>

Summary of lecture 1 - 10.


## Newton's 3 laws {#newton-s-3-laws}

1. You need a force to change a velocity.
Whereas aristoteles would have said, you need a force to move an object.
2. \\( F = ma \\)
This doesn't survive QM.
3. \\( F\_{12} = - F\_{21}\\)


## How we can remove time from the equation. {#how-we-can-remove-time-from-the-equation-dot}

\begin{align\*}
& \frac{dv}{dt} = a \\\\
& v \frac{dv}{dt} = a v \\\\
& \frac{d}{d t}\left(\frac{v^2}{2}\right)=a \frac{d x}{d t} \\\\
& d\left(\frac{v^2}{2}\right)=a d x \\\\
& \int a\left(\frac{v^2}{2}\right)=\int a d x \\\\
& \frac{v^2}{2}-\frac{v\_0^2}{2}=a\left(x-x\_0\right) \\\\
\end{align\*}


## Conservation of energy {#conservation-of-energy}

Energy is work done over distance.

\begin{align\*}
& \frac{v^2}{2}-\frac{v\_0^2}{2}=a d \\\\
& \frac{v^2}{2}-\frac{v\_0^2}{2}=\frac{F}{m} d \\\\
& m\frac{v^2}{2}- m\frac{v\_0^2}{2}=F d = W \text{and we name } m\frac{v^2}{2} \text{ kinetic energy}\\\\
& \sum{\Delta K} = \sum{F \Delta x} \\\\
& K\_{2} - K\_{1} = \int\_{x\_{1}}^{x\_{2}} F(x) dx = G\_{2} - G\_{1} \text{where } \frac{d G}{d x} = F \\\\
& m g h\_{1} - m g h\_{2}
\end{align\*}


## Circular motion {#circular-motion}

We are considering a motion where rotation center is fixed. Then, analogue of distance in linear motion is the angle in the circular motion. After all, radian is a unit of distance (in unit circle).

\begin{alignat\*}{2}
v = \frac{d x}{d t} &\quad& w = \frac{d \theta}{d t} \\\\
p = m v && L = I w \\\\
\end{alignat\*}

\begin{align\*}
& \vec{R}=R(\vec{i} \cos \omega t+\vec{\jmath} \sin \omega t) \\\\
& \frac{d \vec{R}}{d t}=R(-\omega \vec{i} \sin \omega t+ \omega \vec{\jmath} \cos \omega t) \\\\
& \frac{d^{2} \vec{R}}{d t^{2}}=R(-\omega^{2} \vec{i} \cos \omega t - \omega^{2} \vec{\jmath} \sin \omega t) \\\\
& = -\omega^{2} \vec{R}
\end{align\*}


## Torque {#torque}

We have an analogue of \\( F = m a\\) in a circular motion.

\begin{align\*}
& \Delta W = F r \Delta \theta = \Delta K \\\\
& \Delta K = \frac{1}{2} I (W^{2} - W\_{0}^{2}) = \frac{1}{2} I 2 \alpha \Delta \theta \\\\
& F r \Delta \theta = \frac{1}{2} I 2 \alpha \Delta \theta \\\\
& F r = I \alpha
\end{align\*}


## conservation of angular momentum {#conservation-of-angular-momentum}

Without external torque, angular momentum is preserved.

\begin{alignat\*}{2}
p = m v &\quad& L = I w \\\\
\frac{d p}{d t} = 0 &\quad& \frac{d L}{d t} = 0 \\\\
m v\_{1} = m v\_{2} && I w\_{1} = I w\_{2} \\\\
\end{alignat\*}


## Question {#question}

Newton's 2nd law is a conjecture. and Testing it needs to use the conjecture as well.
If we use spring as discussed in the lecture to get the \\(F = -kx\\), deriving it (by experimentation) would also have used the law.

So using the two equation as in the lecture beats the purpose of showing \\( F = ma \\), in a sense that we don't know what the true \\( F \\) is for a string.

\begin{align\*}
F &= ma
F &= -kx
\end{align\*}

What it achieved is that, the \\( F = ma \\) is consistent in describing the motion if we accept the equation as the law.


## Question {#question}

To make an object to orbit, we use \\( ma = m \frac{v^2}{R} \\), so how do we make sure the calculation is exact in a non idealized (with friction) environment?
Do the object has a way of changing its verlocity when needed?


## Question {#question}

in inelastic collision, energy might be lost and conservation of energy might not hold. But momentum is always preserved (when there's no external net force).
Why? well that's what Newton's 1st law is saying. The velocity doesn't change. So 1st law survives.
What happens to internal forces and their effect? they all cancel out.


## Question {#question}

conservation of linear momentum says \\( mv \\) is constant
conservation of angular momentum says \\( Iw \\) is constant
But \\( Iw = m r^2 \* frac{v}{r} = m r v\\), so it says a different quantity is constant.

Now if i were, given an instantaneous information of m and v, I can tell its linear momentum,
but I can't tell the angular momentum unless I need a frame of reference (whree the rotation axis is).

So it feels strange to have two laws where one says \\( mv \\) is constant, and another \\( m r v \\) is constant, and this two statement can't be true at the same time, if the first statement is also conditional statement.

It turns out, you need also need a frame of reference for the linear momentum.
You need an inertial frame of reference for linear momentum.
You need an rotaional(?) frame of reference for angular momentum.

And so it happens there are deeper explanation which I don't yet understand: Noether's theorem. The conservation of linear momentum corresponds to the translational symmetry of space, while the conservation of angular momentum corresponds to its rotational symmetry


## Question {#question}

newton's first law says linear momentum is preserved.

Whereas conservation of angular momentum says the angular momentum is conserved unless there's no external torque.
But then, the rotating object is changing its velocity and according to newton's 2nd law, it needs a force, so what's going on?

Angular motion is another way of saying there's a force acting perpendicular to the motion.
So, angular motion and the related law is describing an object where a there is a perpendicular force acting on it.
And this centripetal force does change the velocity of the object.

So we are considering this special situation where perpendicular force is acting on an object, and linear velocity is changing (either only direction or both direction and magnitude).
But since this force is perpendicular, the perpendicular velocity isn't changing and consequently perpendicular momentum which is angular momentum isnt's changing.

It's easy to see when center of rotation is fixed, can this thought be still regarded true when center of rotation is changing by attributing the change to linear motion so as to keep the r in \\(m r v \\) to remain constant?

on to next question


## Question {#question}

A thought experiment where the mass of the Sun changes but the gravitational force exerted on the earth remains exactly the same due to a corresponding change in the distance between the sun and earth.

If sun changes its mass and distance but keep the gravitational force equal, would earth move to a different path than it would have moved?

Put it differently, does the same force (with same magnitude and same direction) on earth would change the course of Earth's movement when the force is exerted from different distance d1 and d2?

I'd like to guess it doesn't change the path, and that was the core idea of understanding angular momentum and torque, but I can't yet prove it by equation.
(My thought experiment applies same force so earth has to move the same, that's the argument)

If the answer to above question is yes, can we talk about torque without thinking of the axis of rotation? \\( r \\)
No, if there's no axis of rotation, we can't talk about rotation, so torque is not defined.
183 changes: 183 additions & 0 deletions content/posts/2023-10-25-kalman-filter.md
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,183 @@
+++
title = "Kalman Filter"
author = ["littlehome"]
draft = false
+++

## Resources {#resources}

<https://youtu.be/-DiZGpAh7T4?si=YN-a_mZ1gdwlt54n>

<https://web.mit.edu/kirtley/kirtley/binlustuff/literature/control/Kalman%20filter.pdf>

<https://arxiv.org/pdf/1910.03558.pdf>

<https://chat.openai.com/share/d2e080cf-607c-4492-a359-f43c4e7e3daa> (gpt)


## Kalman filter {#kalman-filter}

The Kalman filter operates in a loop with the main stages being prediction and update.

It has model (state) which drives the pheonomenon (observation)
The transition from prev state to next state is linear.
The state to observation is linear.

With updated observation, you update your model so that your prediction is getting better by balancing the weight of 'model prediction' and 'current observation'.

For instance, where you are moving linearly, kalman filter tries to construct a model (parameters) to describe observation better on each observation.

1. ****Initialization****:
- Start with an initial estimate of the state, \\( x\_0 \\).
- Provide an initial estimate of the error covariance, \\( P\_0 \\).

2. ****Prediction Phase****:
- ****State Prediction****: Using the system dynamics, predict the state at the next time step. This prediction uses the most recent state estimate.
\\( x\_{k+1|k} = Ax\_{k|k} \\)
- ****Error Covariance Prediction****: Also predict the error covariance associated with this state prediction.
\\( P\_{k+1|k} = AP\_{k|k}A^T + Q \\)
where \\( Q \\) is the process noise covariance.

3. ****Measurement Update (Correction) Phase****:
- ****Compute Kalman Gain****: This determines how much the prediction should be adjusted based on the new measurement. It's a balance between how much we trust our prediction vs. the measurement.
\\( K\_{k+1} = P\_{k+1|k} H^T (HP\_{k+1|k} H^T + R)^{-1} \\)
where \\( H \\) is the measurement model and \\( R \\) is the measurement noise covariance.
- ****State Estimate Update****: Correct the prediction with the measurement to get a refined state estimate.
\\( x\_{k+1|k+1} = x\_{k+1|k} + K\_{k+1}(y\_{k+1} - Hx\_{k+1|k}) \\)
where \\( y\_{k+1} \\) is the measurement at time \\( k+1 \\).
- ****Error Covariance Update****: Update the error covariance to reflect the new knowledge from the measurement.
\\( P\_{k+1|k+1} = (I - K\_{k+1}H)P\_{k+1|k} \\)

4. ****Loop****:
- Return to the Prediction Phase and repeat for the next time step.

Throughout this process, the Kalman filter effectively combines knowledge from:

- The system's own dynamics (how it believes it evolves over time).
- The new measurements (direct observations of the system).

By combining these two sources of information, the Kalman filter provides a statistically optimal estimate of the state of the system under the assumptions of the model.


### Derivation of Kalman gain {#derivation-of-kalman-gain}

Certainly! The derivation of the Kalman gain can indeed be motivated by a similar logic as the one we followed earlier.

To derive the Kalman gain, we're essentially trying to minimize the posterior error (or the error after the update), given the new measurement. Like before, let's denote the estimation error after incorporating the measurement as:

\\[ e\_{k+1|k+1} = x\_{k+1} - x\_{k+1|k+1} \\]

With the Kalman update, the updated estimate is:

\\[ x\_{k+1|k+1} = x\_{k+1|k} + K\_{k+1}(y\_{k+1} - C x\_{k+1|k}) \\]

Substituting in for our error:

\\[ e\_{k+1|k+1} = x\_{k+1} - x\_{k+1|k} - K\_{k+1}(y\_{k+1} - C x\_{k+1|k}) \\]

Now, using our system dynamics, we know:
\\[ y\_{k+1} = C x\_{k+1} + v\_{k+1} \\]
Where \\( v\_{k+1} \\) is our measurement noise.

Rearranging gives:
\\[ e\_{k+1|k+1} = (I - K\_{k+1} C) e\_{k+1|k} - K\_{k+1} v\_{k+1} \\]

To minimize the expected value of the squared error, we can consider its covariance matrix:

\\[ E[e\_{k+1|k+1} e\_{k+1|k+1}^T] \\]

Expanding this out and considering the cross term of the error and measurement noise, we get:
\\[ E[e\_{k+1|k} e\_{k+1|k}^T] - E[K\_{k+1} C e\_{k+1|k} e\_{k+1|k}^T] - E[e\_{k+1|k} e\_{k+1|k}^T C^T K\_{k+1}^T] + E[K\_{k+1} v\_{k+1} v\_{k+1}^T K\_{k+1}^T] \\]

We aim to minimize this expression by selecting the best \\( K\_{k+1} \\). Differentiating w.r.t \\( K\_{k+1} \\) and setting the result to zero (like we did previously), and recognizing the error covariance matrix \\( P\_{k+1|k} = E[e\_{k+1|k} e\_{k+1|k}^T] \\), we get:

\\[ K\_{k+1} = P\_{k+1|k} C^T (C P\_{k+1|k} C^T + R)^{-1} \\]

This is the Kalman gain. The process we followed mirrors the reasoning we used in the simpler estimation problem, but now in the context of state-space models and dynamical systems.


## Implementaion {#implementaion}

```python
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(0)
```

```python
def get_volt():
"""Measure voltage."""
v = np.random.normal(0, 2) # v: measurement noise.
volt_true = 14.4 # volt_true: True voltage [V].
measured = volt_true + v # z_volt_meas: Measured Voltage [V] (observable).
return measured


def kalman_filter(measured, u_estimate, P):
"""Kalman Filter Algorithm for One Variable."""
# (1) Prediction.
# u_estimate is previous step's estimate
# u_prediction is this step's prediction
# P_prediction is covariance (?)

u_prediction = A * u_estimate # prediction by last estimate (best guess)
P_prediction = A * P * A + Q # update uncertainties

# (2) Kalman Gain.
K = P_prediction * H / (H * P_prediction * H + R)

# (3) Estimation.
u_estimate = u_prediction + K * (measured - H * u_prediction)

# (4) Error Covariance.
P = P_prediction - K * H * P_prediction

return u_estimate, P
```

```python
# Input parameters.
time_end = 10
dt = 0.2

# Initialization for system model.
A = 1
H = 1
Q = 0
R = 4
# Initialization for estimation.
x_0 = 12 # 14 for book.
P_0 = 6

time = np.arange(0, time_end, dt)
n_samples = len(time)
volt_measured_save = np.zeros(n_samples)
volt_esti_save = np.zeros(n_samples)
```

```python
u_estimate, P = None, None
for i in range(n_samples):
measured = get_volt()
if i == 0:
u_estimate, P = x_0, P_0
else:
u_estimate, P = kalman_filter(measured, u_estimate, P)

volt_measured_save[i] = measured
volt_esti_save[i] = u_estimate
```

```python
plt.plot(time, volt_measured_save, 'r*--', label='Measurements')
plt.plot(time, volt_esti_save, 'bo-', label='Kalman Filter')
plt.legend(loc='upper left')
plt.title('Measurements v.s. Estimation (Kalman Filter)')
plt.xlabel('Time [sec]')
plt.ylabel('Voltage [V]')
# plt.savefig('png/simple_kalman_filter.png')
plt.show()
```

{{< figure src="/ox-hugo/2023-10-25-kalman.png" >}}
Binary file added static/ox-hugo/2023-10-25-kalman.png
View file
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.

0 comments on commit 78a98a9

Please sign in to comment.