```python
# Boilerplate imports
from numpy import array as ary; from numpy import log as ln
from numpy import cos, sin, pi, sqrt, exp, arccos, arcsin
tau = 2*pi
import numpy as np;
from matplotlib import pyplot as plt
# Linear algebra functions
from numpy.linalg import inv, pinv, det, eig, eigh, eigvals
from matplotlib.patches import Ellipse # for plotting ellipse
from collections import namedtuple
```


```python
# program to generate a covariance matrix where the variance values are fixed at [2,2].
def generate_cov(off_diag):
    cov = ary([[2.0, off_diag], [off_diag, 2.0]])
    return cov

Dots = namedtuple('Dots', ['points', 'area'])
def get_encircled_dots(covariance_matrix, bounds=[-10, 10], resolution=300):
    pt_list = np.linspace(*bounds, resolution)
    unit_square_size = ((bounds[1]-bounds[0])/resolution)**2
    points = ary(np.meshgrid(pt_list, pt_list)).T.flatten().reshape([-1,2])
    
    chi2_level = ((points @ inv(covariance_matrix)) * points).sum(axis=1)
    mask = chi2_level <= 1 # choose only points within the error ellipse
    area = sum(mask)*unit_square_size
    
    return Dots(points[mask], area)


PLOT_CIRCLE = True
# if PLOT_CIRCLE: plot the error ellipse using 11 differernt covariance values;
# else: plot plot the variation of the error ellipse area wrt. the covariance value.
if not PLOT_CIRCLE:
    determinant_cov, determinant_inv_cov, size = [], [], []
```

Setup:

- We have a point cloud with a center at (0,0).
- 1 sigma of the points (68% of them) lies within the error ellipse. (We won't be plotting the remaining 32%)
- The top left element of the covariance matrix describes the width of this ellipse, and the bottom right describes the height of the ellipse.
- Therefore, varying the covariance value (i.e. the symmetric off-diagonal terms) should only make the ellipse into a thinner ellipse that is leaning left/right.


```python
# plot the error ellipse
cov_range = np.linspace(-1.98, 1.98, 11)
for i in cov_range:
    cov = generate_cov(i)
    (minor, major) , eig_mat = eig(inv(cov) * det(cov))
    mat = inv(eig_mat)
    orientation = np.mean([arcsin(mat[0,1]), -arcsin(mat[1,0])])*np.sign(mat[0,0])
    fig, ax = plt.subplots()
    ellipse = Ellipse([0,0], # centered at the origin
                2*sqrt(major),
                2*sqrt(minor),
                np.rad2deg(orientation)
                )
    # DUUUUDE I got the width=2*sqrt(major)/sqrt(det(inv(cov))) equation by trial and error LMAO
    ax.add_patch(ellipse)
    ax.scatter(*(get_encircled_dots(cov).points.T), marker='+', alpha=0.4, color='C1', zorder=10) # scatter plot approach
    plt.show()

```

Using the code block above we can verify that we have correctly plotted the error ellipse correctly: All points with $\chi^2 \le 1$ are plotted; and it concides exactly with the error ellipse.
The parameters about the ellipse is closely related to the following matrix:

\begin{equation}
M = S^{-1} \dot det(S)
\end{equation}

where S is the covariance matrix

The major radius is equal to sqrt(the larger eigen value of $M$), minor radius is equal to sqrt(the smaller eigen value of $M$).
To draw the ellipse on our graph, we first draw an ellipse with those specified radii (major axis in the horizontal direction), then apply the rotation matrix as described by the eigenvector matrix of $M$.

We can also show that the covariance ellipse area is equal to the expression $det(S)$ (i.e. determinant of the covariance matrix).

The ellipse area was empirically calculated by counting the number of points that the ellipse covers when spread over an evenly spaced grid.


```python
determinant_cov, size = [], []
for i in np.linspace(-0, 1.98, 30):
    cov = generate_cov(i)
    determinant_cov.append(det(cov))
    size.append(get_encircled_dots(cov).area)

plt.plot(cov_range, sqrt(determinant_cov)*pi, label='cov determinant')
plt.xlabel('covariance (off-diagonal elements) value')
plt.ylabel('area/area prediction/other quantities')
plt.plot(cov_range, size, label='ellipse size')
plt.legend()
plt.show()

print("Notice that these two lines overlap very well.")
print("In fact, they would be exactly the same if the number of samples we take approaches infinity.")
```

The important conclusion from this project is that, if we fix the variance values, increasing the absolute value of the covariance will make the error ellipse thinner. And the specific algorithm required to plot the error ellipse is also found.

17005e240b8b2511a21eaf6b5e067c9ba8f2f153
ipynb
Jupyter Notebook
2d.ipynb
OceanNuclear/Covariance
9b3c8ed4cf01fe36824de26b414ce7386d1f87f5
MIT
2d.ipynb
OceanNuclear/Covariance
9b3c8ed4cf01fe36824de26b414ce7386d1f87f5
MIT
2d.ipynb
OceanNuclear/Covariance
9b3c8ed4cf01fe36824de26b414ce7386d1f87f5
MIT
Qwen/Qwen-72B
1. YES
2. YES