adding more material, miss theorem and its proof
This commit is contained in:
mhjensen
+136
-16
@@ -341,7 +341,7 @@ We have a data set defined by a design/feature matrix $\bm{X}$ (see below for it
|
||||
===== Introducing the Covariance and Correlation functions =====
|
||||
|
||||
Suppose we have defined two vectors
|
||||
$\hat{x} and \hat{y} with $n$ elements each. The covariance matrix $\bm{C}is defined as
|
||||
$\hat{x}$ and $\hat{y}$ with $n$ elements each. The covariance matrix $\bm{C}$ is defined as
|
||||
!bt
|
||||
\[
|
||||
\bm{C}[\bm{x},\bm{y}] = \begin{bmatrix} cov[\bm{x},\bm{x}] & cov[\bm{x},\bm{y}] \\
|
||||
@@ -378,7 +378,7 @@ correlation function
|
||||
|
||||
!bt
|
||||
\[
|
||||
corr[\bm{x},\bm{y}]=\frac{cov[\bm{x},\bm{y}]}{\sqrt{var[\bm{x}]\var\bm{y}]}}.
|
||||
corr[\bm{x},\bm{y}]=\frac{cov[\bm{x},\bm{y}]}{\sqrt{var[\bm{x}]\var[\bm{y}]}}.
|
||||
\]
|
||||
!et
|
||||
|
||||
@@ -456,14 +456,20 @@ corr[\bm{x}_{p-1},\bm{x}_0] & corr[\bm{x}_{p-1},\bm{x}_1] & corr[\bm{x}_{p-1},
|
||||
!et
|
||||
|
||||
|
||||
!split
|
||||
===== Covariance Matrix Examples =====
|
||||
|
||||
|
||||
The Numpy function _np.cov_ calculates the covariance elements using
|
||||
the factor $1/(n-1)$ instead of $1/n$ since it assumes we do not have
|
||||
the exact mean values. The following simple function uses the
|
||||
_np.vstack_ function which takes each vector of dimension $1\times n$
|
||||
and produces a $2\times n$ matrix $\bm{W}$
|
||||
|
||||
|
||||
The Numpy function _np.cov_ calculates the covariance elements using the factor $1/(n-1)$ instead of $1/n$ since it assumes we do not have the exact mean values.
|
||||
The following simple function uses the _np.vstack_ function which takes each vector of dimension $1\times n$ and produces a $2\times n$ matrix $\hat{W}$
|
||||
!bt
|
||||
\[
|
||||
\hat{W} = \begin{bmatrix} x_0 & y_0 \\
|
||||
\bm{W} = \begin{bmatrix} x_0 & y_0 \\
|
||||
x_1 & y_1 \\
|
||||
x_2 & y_2\\
|
||||
\dots & \dots \\
|
||||
@@ -473,30 +479,144 @@ The following simple function uses the _np.vstack_ function which takes each vec
|
||||
\]
|
||||
!et
|
||||
|
||||
which in turn is converted into into the $3\times 3$ covariance matrix
|
||||
$\hat{\Sigma}$ via the Numpy function _np.cov()_. We note that we can also calculate
|
||||
the mean value of each set of samples $\hat{x}$ etc using the Numpy
|
||||
which in turn is converted into into the $2\times 2$ covariance matrix
|
||||
$\bm{C}$ via the Numpy function _np.cov()_. We note that we can also calculate
|
||||
the mean value of each set of samples $\bm{x}$ etc using the Numpy
|
||||
function _np.mean(x)_. We can also extract the eigenvalues of the
|
||||
covariance matrix through the _np.linalg.eig()_ function.
|
||||
|
||||
!bc pycod
|
||||
# Importing various packages
|
||||
import numpy as np
|
||||
|
||||
n = 100
|
||||
x = np.random.normal(size=n)
|
||||
print(np.mean(x))
|
||||
y = 4+3*x+np.random.normal(size=n)
|
||||
print(np.mean(y))
|
||||
z = x**3+np.random.normal(size=n)
|
||||
print(np.mean(z))
|
||||
W = np.vstack((x, y, z))
|
||||
Sigma = np.cov(W)
|
||||
print(Sigma)
|
||||
Eigvals, Eigvecs = np.linalg.eig(Sigma)
|
||||
print(Eigvals)
|
||||
W = np.vstack((x, y))
|
||||
C = np.cov(W)
|
||||
print(C)
|
||||
!ec
|
||||
|
||||
!split
|
||||
===== Correlation Matrix =====
|
||||
|
||||
The previous example can be converted into the correlation matrix by
|
||||
simply scaling the matrix elements with the variances. We should also
|
||||
subtract the mean values for each column. This leads to the following
|
||||
code which sets up the correlations matrix for the previous example in
|
||||
a more brute force way. Here we scale the mean values for each column of the design matrix, calculate the relevant mean values and variances and then finally set up the $2\times 2$ correlation matrix (since we have only two vectors).
|
||||
|
||||
!bc pycod
|
||||
import numpy as np
|
||||
n = 100
|
||||
# define two vectors
|
||||
x = np.random.random(size=n)
|
||||
y = 4+3*x+np.random.normal(size=n)
|
||||
#scaling the x and y vectors
|
||||
x = x - np.mean(x)
|
||||
y = y - np.mean(y)
|
||||
variance_x = np.sum(x@x)/n
|
||||
variance_y = np.sum(y@y)/n
|
||||
print(variance_x)
|
||||
print(variance_y)
|
||||
cov_xy = np.sum(x@y)/n
|
||||
cov_xx = np.sum(x@x)/n
|
||||
cov_yy = np.sum(y@y)/n
|
||||
C = np.zeros((2,2))
|
||||
C[0,0]= cov_xx/variance_x
|
||||
C[1,1]= cov_yy/variance_y
|
||||
C[0,1]= cov_xy/np.sqrt(variance_y*variance_x)
|
||||
C[1,0]= C[0,1]
|
||||
print(C)
|
||||
!ec
|
||||
|
||||
We see that the matrix elements along the diagonal are one as they
|
||||
should be and that the matrix is symmetric. Furthermore, diagonalizing
|
||||
this matrix we easily see that it is a positive definite matrix.
|
||||
|
||||
The above procedure with _numpy_ can be made more compact if we use _pandas_.
|
||||
|
||||
!split
|
||||
===== Correlation Matrix with Pandas =====
|
||||
|
||||
We whow here how we can set up the correlation matrix using _pandas_, as done in this simple code
|
||||
!bc pycod
|
||||
import numpy as np
|
||||
import pandas as pd
|
||||
n = 10
|
||||
x = np.random.normal(size=n)
|
||||
x = x - np.mean(x)
|
||||
y = 4+3*x+np.random.normal(size=n)
|
||||
y = y - np.mean(y)
|
||||
X = (np.vstack((x, y))).T
|
||||
print(X)
|
||||
Xpd = pd.DataFrame(X)
|
||||
print(Xpd)
|
||||
correlation_matrix = Xpd.corr()
|
||||
print(correlation_matrix)
|
||||
!ec
|
||||
|
||||
|
||||
We expand this model to the Franke function discussed above.
|
||||
|
||||
!split
|
||||
===== Correlation Matrix with Pandas and the Franke function =====
|
||||
|
||||
!bc pycod
|
||||
# Common imports
|
||||
import numpy as np
|
||||
import pandas as pd
|
||||
|
||||
|
||||
def FrankeFunction(x,y):
|
||||
term1 = 0.75*np.exp(-(0.25*(9*x-2)**2) - 0.25*((9*y-2)**2))
|
||||
term2 = 0.75*np.exp(-((9*x+1)**2)/49.0 - 0.1*(9*y+1))
|
||||
term3 = 0.5*np.exp(-(9*x-7)**2/4.0 - 0.25*((9*y-3)**2))
|
||||
term4 = -0.2*np.exp(-(9*x-4)**2 - (9*y-7)**2)
|
||||
return term1 + term2 + term3 + term4
|
||||
|
||||
|
||||
def create_X(x, y, n ):
|
||||
if len(x.shape) > 1:
|
||||
x = np.ravel(x)
|
||||
y = np.ravel(y)
|
||||
|
||||
N = len(x)
|
||||
l = int((n+1)*(n+2)/2) # Number of elements in beta
|
||||
X = np.ones((N,l))
|
||||
|
||||
for i in range(1,n+1):
|
||||
q = int((i)*(i+1)/2)
|
||||
for k in range(i+1):
|
||||
X[:,q+k] = (x**(i-k))*(y**k)
|
||||
|
||||
return X
|
||||
|
||||
|
||||
# Making meshgrid of datapoints and compute Franke's function
|
||||
n = 4
|
||||
N = 100
|
||||
x = np.sort(np.random.uniform(0, 1, N))
|
||||
y = np.sort(np.random.uniform(0, 1, N))
|
||||
z = FrankeFunction(x, y)
|
||||
X = create_X(x, y, n=n)
|
||||
|
||||
Xpd = pd.DataFrame(X)
|
||||
# subtract the mean values and set up the covariance matrix
|
||||
Xpd = Xpd - Xpd.mean()
|
||||
covariance_matrix = Xpd.cov()
|
||||
print(covariance_matrix)
|
||||
!ec
|
||||
|
||||
We note here that the covariance is zero for the first rows and
|
||||
columns since all matrix elements in the design matrix were set to one
|
||||
(we are fitting the function in terms of a polynomial of degree $n$).
|
||||
|
||||
This means that the variance for these elements will be zero and will
|
||||
cause problems when we set up the correlation matrix. We can simply
|
||||
drop these elements as follows and then construct the correlation
|
||||
matrix.
|
||||
|
||||
|
||||
|
||||
|
||||
@@ -0,0 +1,24 @@
|
||||
# Importing various packages
|
||||
import numpy as np
|
||||
n = 100
|
||||
# define two vectors
|
||||
x = np.random.random(size=n)
|
||||
y = 4+3*x+np.random.normal(size=n)
|
||||
#scaling the x and y vectors
|
||||
x = x - np.mean(x)
|
||||
y = y - np.mean(y)
|
||||
variance_x = np.sum(x@x)/n
|
||||
variance_y = np.sum(y@y)/n
|
||||
print(variance_x)
|
||||
print(variance_y)
|
||||
cov_xy = np.sum(x@y)/n
|
||||
cov_xx = np.sum(x@x)/n
|
||||
cov_yy = np.sum(y@y)/n
|
||||
C = np.zeros((2,2))
|
||||
C[0,0]= cov_xx/variance_x
|
||||
C[1,1]= cov_yy/variance_y
|
||||
C[0,1]= cov_xy/np.sqrt(variance_y*variance_x)
|
||||
C[1,0]= C[0,1]
|
||||
print(C)
|
||||
Eigvals, Eigvecs = np.linalg.eig(C)
|
||||
print(Eigvals)
|
||||
@@ -0,0 +1,43 @@
|
||||
# Common imports
|
||||
import numpy as np
|
||||
import pandas as pd
|
||||
|
||||
|
||||
def FrankeFunction(x,y):
|
||||
term1 = 0.75*np.exp(-(0.25*(9*x-2)**2) - 0.25*((9*y-2)**2))
|
||||
term2 = 0.75*np.exp(-((9*x+1)**2)/49.0 - 0.1*(9*y+1))
|
||||
term3 = 0.5*np.exp(-(9*x-7)**2/4.0 - 0.25*((9*y-3)**2))
|
||||
term4 = -0.2*np.exp(-(9*x-4)**2 - (9*y-7)**2)
|
||||
return term1 + term2 + term3 + term4
|
||||
|
||||
|
||||
def create_X(x, y, n ):
|
||||
if len(x.shape) > 1:
|
||||
x = np.ravel(x)
|
||||
y = np.ravel(y)
|
||||
|
||||
N = len(x)
|
||||
l = int((n+1)*(n+2)/2) # Number of elements in beta
|
||||
X = np.ones((N,l))
|
||||
|
||||
for i in range(1,n+1):
|
||||
q = int((i)*(i+1)/2)
|
||||
for k in range(i+1):
|
||||
X[:,q+k] = (x**(i-k))*(y**k)
|
||||
|
||||
return X
|
||||
|
||||
|
||||
# Making meshgrid of datapoints and compute Franke's function
|
||||
n = 4
|
||||
N = 1000
|
||||
x = np.sort(np.random.uniform(0, 1, N))
|
||||
y = np.sort(np.random.uniform(0, 1, N))
|
||||
z = FrankeFunction(x, y)
|
||||
X = create_X(x, y, n=n)
|
||||
|
||||
Xpd = pd.DataFrame(X)
|
||||
Xpd = Xpd - Xpd.mean()
|
||||
correlation_matrix = Xpd.cov()
|
||||
print(correlation_matrix)
|
||||
|
||||
@@ -0,0 +1,28 @@
|
||||
# Importing various packages
|
||||
import numpy as np
|
||||
n = 10
|
||||
x = np.random.normal(size=n)
|
||||
x = x - np.mean(x)
|
||||
y = 4+3*x+np.random.normal(size=n)
|
||||
y = y - np.mean(y)
|
||||
X = (np.vstack((x, y))).T
|
||||
print(X)
|
||||
import pandas as pd
|
||||
Xpd = pd.DataFrame(X)
|
||||
print(Xpd)
|
||||
correlation_matrix = Xpd.corr()
|
||||
print(correlation_matrix)
|
||||
|
||||
|
||||
|
||||
variance_x = np.sum(x@x)/n
|
||||
variance_y = np.sum(y@y)/n
|
||||
cov_xy = np.sum(x@y)/n
|
||||
cov_xx = np.sum(x@x)/n
|
||||
cov_yy = np.sum(y@y)/n
|
||||
C = np.zeros((2,2))
|
||||
C[0,0]= cov_xx/variance_x
|
||||
C[1,1]= cov_yy/variance_y
|
||||
C[0,1]= cov_xy/np.sqrt(variance_y*variance_x)
|
||||
C[1,0]= C[0,1]
|
||||
print(C)
|
||||
Reference in New Issue
Block a user