typos in reg slides

This commit is contained in:
mhjensen
2020-09-11 11:10:01 +02:00
parent bea41d2818
commit 243a68c942
2 changed files with 252 additions and 225 deletions
+201 -176
View File
@@ -381,9 +381,7 @@
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline\n",
@@ -888,9 +886,7 @@
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# matrix inversion to find beta\n",
@@ -909,9 +905,7 @@
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"fit = np.linalg.lstsq(X, Energies, rcond =None)[0]\n",
@@ -928,9 +922,7 @@
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"Masses['Eapprox'] = ytilde\n",
@@ -960,9 +952,7 @@
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"def R2(y_data, y_model):\n",
@@ -979,9 +969,7 @@
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"print(R2(Energies,ytilde))"
@@ -997,9 +985,7 @@
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"def MSE(y_data,y_model):\n",
@@ -1019,9 +1005,7 @@
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"def RelativeError(y_data,y_model):\n",
@@ -1398,9 +1382,7 @@
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Common imports\n",
@@ -1519,9 +1501,7 @@
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import os\n",
@@ -1638,9 +1618,7 @@
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -1660,9 +1638,7 @@
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from sklearn.datasets import load_boston\n",
@@ -1684,9 +1660,7 @@
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"boston = pd.DataFrame(boston_dataset.data, columns=boston_dataset.feature_names)\n",
@@ -1704,9 +1678,7 @@
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# check for missing values in all the columns\n",
@@ -1723,9 +1695,7 @@
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# set the size of the figure\n",
@@ -1746,9 +1716,7 @@
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# compute the pair wise correlation for all columns \n",
@@ -1768,9 +1736,7 @@
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"plt.figure(figsize=(20, 5))\n",
@@ -1798,9 +1764,7 @@
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"X = pd.DataFrame(np.c_[boston['LSTAT'], boston['RM']], columns = ['LSTAT','RM'])\n",
@@ -1817,9 +1781,7 @@
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from sklearn.model_selection import train_test_split\n",
@@ -1843,9 +1805,7 @@
{
"cell_type": "code",
"execution_count": 20,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from sklearn.linear_model import LinearRegression\n",
@@ -1884,9 +1844,7 @@
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# plotting the y_test vs y_pred\n",
@@ -1971,9 +1929,7 @@
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Common imports\n",
@@ -2839,11 +2795,33 @@
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"collapsed": false
},
"outputs": [],
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[[ 1. -1. 2.]\n",
" [ 1. 0. 1.]\n",
" [ 1. 2. -1.]\n",
" [ 1. 1. 0.]]\n",
"[[ 4. 2. 2.]\n",
" [ 2. 6. -4.]\n",
" [ 2. -4. 6.]]\n",
"[[-9.57425734e-17 8.16496581e-01 -5.77350269e-01]\n",
" [-7.07106781e-01 4.08248290e-01 5.77350269e-01]\n",
" [ 7.07106781e-01 4.08248290e-01 5.77350269e-01]]\n",
"[1.00000000e+01 6.00000000e+00 9.10898112e-32]\n",
"[[ 3.33066907e-17 -7.07106781e-01 7.07106781e-01]\n",
" [ 8.16496581e-01 4.08248290e-01 4.08248290e-01]\n",
" [ 5.77350269e-01 -5.77350269e-01 -5.77350269e-01]]\n",
"[[-3.65939208e+30 3.65939208e+30 3.65939208e+30]\n",
" [ 3.65939208e+30 -3.65939208e+30 -3.65939208e+30]\n",
" [ 3.65939208e+30 -3.65939208e+30 -3.65939208e+30]]\n"
]
}
],
"source": [
"import numpy as np\n",
"# SVD inversion\n",
@@ -3198,11 +3176,20 @@
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"collapsed": false
},
"outputs": [],
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.06921280734901214\n",
"4.282897881229527\n",
"[[ 1.16689268 3.45320489]\n",
" [ 3.45320489 11.11096825]]\n"
]
}
],
"source": [
"# Importing various packages\n",
"import numpy as np\n",
@@ -3231,11 +3218,20 @@
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {
"collapsed": false
},
"outputs": [],
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.07629059232706191\n",
"1.3689323216991471\n",
"[[1. 0.63967543]\n",
" [0.63967543 1. ]]\n"
]
}
],
"source": [
"import numpy as np\n",
"n = 100\n",
@@ -3277,11 +3273,40 @@
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {
"collapsed": false
},
"outputs": [],
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[[-0.55783318 -2.36855737]\n",
" [-0.48038212 -1.12123656]\n",
" [-0.90619932 -1.43972046]\n",
" [ 0.0509597 0.25649197]\n",
" [-0.00589746 1.0155279 ]\n",
" [ 0.42557278 1.62052588]\n",
" [ 0.76914274 -1.09850401]\n",
" [-0.13134882 0.66162568]\n",
" [ 0.85273632 2.51549691]\n",
" [-0.01675064 -0.04164992]]\n",
" 0 1\n",
"0 -0.557833 -2.368557\n",
"1 -0.480382 -1.121237\n",
"2 -0.906199 -1.439720\n",
"3 0.050960 0.256492\n",
"4 -0.005897 1.015528\n",
"5 0.425573 1.620526\n",
"6 0.769143 -1.098504\n",
"7 -0.131349 0.661626\n",
"8 0.852736 2.515497\n",
"9 -0.016751 -0.041650\n",
" 0 1\n",
"0 1.000000 0.657956\n",
"1 0.657956 1.000000\n"
]
}
],
"source": [
"import numpy as np\n",
"import pandas as pd\n",
@@ -3309,11 +3334,49 @@
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {
"collapsed": false
},
"outputs": [],
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" 0 1 2 3 4 5 6 7 \\\n",
"0 0.0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 \n",
"1 0.0 0.098274 0.088019 0.102816 0.095130 0.088038 0.095734 0.089350 \n",
"2 0.0 0.088019 0.080153 0.091357 0.085370 0.079873 0.085008 0.079956 \n",
"3 0.0 0.102816 0.091357 0.113039 0.104455 0.096480 0.108457 0.101346 \n",
"4 0.0 0.095130 0.085370 0.104455 0.097137 0.090339 0.100399 0.094305 \n",
"5 0.0 0.088038 0.079873 0.096480 0.090339 0.084640 0.092879 0.087723 \n",
"6 0.0 0.095734 0.085008 0.108457 0.100399 0.092879 0.106219 0.099527 \n",
"7 0.0 0.089350 0.079956 0.101346 0.094305 0.087723 0.099527 0.093665 \n",
"8 0.0 0.083518 0.075349 0.094820 0.088708 0.082989 0.093363 0.088259 \n",
"9 0.0 0.078178 0.071145 0.088809 0.083553 0.078629 0.087665 0.083254 \n",
"10 0.0 0.087401 0.077865 0.101048 0.093860 0.087128 0.100476 0.094478 \n",
"11 0.0 0.081937 0.073490 0.094954 0.088611 0.082657 0.094714 0.089415 \n",
"12 0.0 0.076958 0.069504 0.089383 0.083807 0.078560 0.089433 0.084768 \n",
"13 0.0 0.072416 0.065869 0.084280 0.079402 0.074802 0.084581 0.080492 \n",
"14 0.0 0.068264 0.062550 0.079595 0.075355 0.071347 0.080112 0.076548 \n",
"\n",
" 8 9 10 11 12 13 14 \n",
"0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 \n",
"1 0.083518 0.078178 0.087401 0.081937 0.076958 0.072416 0.068264 \n",
"2 0.075349 0.071145 0.077865 0.073490 0.069504 0.065869 0.062550 \n",
"3 0.094820 0.088809 0.101048 0.094954 0.089383 0.084280 0.079595 \n",
"4 0.088708 0.083553 0.093860 0.088611 0.083807 0.079402 0.075355 \n",
"5 0.082989 0.078629 0.087128 0.082657 0.078560 0.074802 0.071347 \n",
"6 0.093363 0.087665 0.100476 0.094714 0.089433 0.084581 0.080112 \n",
"7 0.088259 0.083254 0.094478 0.089415 0.084768 0.080492 0.076548 \n",
"8 0.083545 0.079175 0.088938 0.084512 0.080443 0.076693 0.073228 \n",
"9 0.079175 0.075391 0.083801 0.079956 0.076416 0.073148 0.070125 \n",
"10 0.088938 0.083801 0.096178 0.090983 0.086211 0.081816 0.077756 \n",
"11 0.084512 0.079956 0.090983 0.086383 0.082149 0.078242 0.074627 \n",
"12 0.080443 0.076416 0.086211 0.082149 0.078403 0.074940 0.071730 \n",
"13 0.076693 0.073148 0.081816 0.078242 0.074940 0.071881 0.069039 \n",
"14 0.073228 0.070125 0.077756 0.074627 0.071730 0.069039 0.066534 \n"
]
}
],
"source": [
"# Common imports\n",
"import numpy as np\n",
@@ -3956,9 +4019,7 @@
{
"cell_type": "code",
"execution_count": 28,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Importing various packages\n",
@@ -3981,9 +4042,7 @@
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -5146,9 +5205,7 @@
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from numpy import *\n",
@@ -5289,9 +5346,7 @@
{
"cell_type": "code",
"execution_count": 31,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"from numpy import *\n",
@@ -5437,9 +5492,7 @@
{
"cell_type": "code",
"execution_count": 32,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -5668,9 +5721,7 @@
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
@@ -5739,9 +5790,7 @@
{
"cell_type": "code",
"execution_count": 34,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
@@ -5840,9 +5889,7 @@
{
"cell_type": "code",
"execution_count": 35,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"\"\"\"\n",
@@ -5928,9 +5975,7 @@
{
"cell_type": "code",
"execution_count": 36,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Common imports\n",
@@ -6024,9 +6069,7 @@
{
"cell_type": "code",
"execution_count": 37,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"# Common imports\n",
@@ -6108,9 +6151,7 @@
{
"cell_type": "code",
"execution_count": 38,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -6191,9 +6232,7 @@
{
"cell_type": "code",
"execution_count": 39,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -6306,9 +6345,7 @@
{
"cell_type": "code",
"execution_count": 40,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"X = np.zeros((n, L ** 2))\n",
@@ -6372,9 +6409,7 @@
{
"cell_type": "code",
"execution_count": 41,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"X_train_own = np.concatenate(\n",
@@ -6390,9 +6425,7 @@
{
"cell_type": "code",
"execution_count": 42,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"def ols_inv(x: np.ndarray, y: np.ndarray) -> np.ndarray:\n",
@@ -6477,9 +6510,7 @@
{
"cell_type": "code",
"execution_count": 43,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"def ols_svd(x: np.ndarray, y: np.ndarray) -> np.ndarray:\n",
@@ -6490,9 +6521,7 @@
{
"cell_type": "code",
"execution_count": 44,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"beta = ols_svd(X_train_own,y_train)"
@@ -6508,9 +6537,7 @@
{
"cell_type": "code",
"execution_count": 45,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"J = beta[1:].reshape(L, L)"
@@ -6526,9 +6553,7 @@
{
"cell_type": "code",
"execution_count": 46,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"fig = plt.figure(figsize=(20, 14))\n",
@@ -6594,9 +6619,7 @@
{
"cell_type": "code",
"execution_count": 47,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
@@ -6702,9 +6725,7 @@
{
"cell_type": "code",
"execution_count": 48,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"X = np.zeros((n, L ** 2))\n",
@@ -6734,9 +6755,7 @@
{
"cell_type": "code",
"execution_count": 49,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"clf = skl.LinearRegression().fit(X_train, y_train)"
@@ -6752,9 +6771,7 @@
{
"cell_type": "code",
"execution_count": 50,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"J_sk = clf.coef_.reshape(L, L)"
@@ -6770,9 +6787,7 @@
{
"cell_type": "code",
"execution_count": 51,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"fig = plt.figure(figsize=(20, 14))\n",
@@ -6828,9 +6843,7 @@
{
"cell_type": "code",
"execution_count": 52,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"_lambda = 0.1\n",
@@ -6881,9 +6894,7 @@
{
"cell_type": "code",
"execution_count": 53,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"clf_lasso = skl.Lasso(alpha=_lambda).fit(X_train, y_train)\n",
@@ -6917,9 +6928,7 @@
{
"cell_type": "code",
"execution_count": 54,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"lambdas = np.logspace(-4, 5, 10)\n",
@@ -6982,9 +6991,7 @@
{
"cell_type": "code",
"execution_count": 55,
"metadata": {
"collapsed": false
},
"metadata": {},
"outputs": [],
"source": [
"fig = plt.figure(figsize=(20, 14))\n",
@@ -7029,7 +7036,25 @@
]
}
],
"metadata": {},
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
+51 -49
View File
@@ -1312,7 +1312,7 @@ This may
however not the be case in general and a standard matrix inversion
algorithm based on say LU, QR or Cholesky decomposition may lead to singularities. We will see examples of this below.
There is however a way to partially circumvent this problem and also gain some insight about the ordinary least squares approach.
There is however a way to partially circumvent this problem and also gain some insights about the ordinary least squares approach, and later shrinkage methods like Ridge and Lasso regressions.
This is given by the _Singular Value Decomposition_ algorithm, perhaps
the most powerful linear algebra algorithm. Let us look at a
@@ -1519,6 +1519,56 @@ If $p > n$, then only the first $n$ columns of $\bm{V}$ are computed and $\bm{\S
The $n=p$ case is obvious, we retain the full SVD.
In general the economy-size SVD leads to less FLOPS and still conserving the desired accuracy.
!split
===== Codes for the SVD =====
!bc pycod
import numpy as np
# SVD inversion
def SVDinv(A):
''' Takes as input a numpy matrix A and returns inv(A) based on singular value decomposition (SVD).
SVD is numerically more stable than the inversion algorithms provided by
numpy and scipy.linalg at the cost of being slower.
'''
U, s, VT = np.linalg.svd(A)
# print('test U')
# print( (np.transpose(U) @ U - U @np.transpose(U)))
# print('test VT')
# print( (np.transpose(VT) @ VT - VT @np.transpose(VT)))
print(U)
print(s)
print(VT)
D = np.zeros((len(U),len(VT)))
for i in range(0,len(VT)):
D[i,i]=s[i]
UT = np.transpose(U); V = np.transpose(VT); invD = np.linalg.inv(D)
return np.matmul(V,np.matmul(invD,UT))
X = np.array([ [1.0, -1.0, 2.0], [1.0, 0.0, 1.0], [1.0, 2.0, -1.0], [1.0, 1.0, 0.0] ])
print(X)
A = np.transpose(X) @ X
print(A)
# Brute force inversion of super-collinear matrix
#B = np.linalg.inv(A)
#print(B)
C = SVDinv(A)
print(C)
!ec
The matrix $\bm{X}$ has columns that are linearly dependent. The first
column is the row-wise sum of the other two columns. The rank of a
matrix (the column rank) is the dimension of space spanned by the
column vectors. The rank of the matrix is the number of linearly
independent columns, in this case just $2$. We see this from the
singular values when running the above code. Running the standard
inversion algorithm for matrix inversion with $\bm{X}^T\bm{X}$ results
in the program terminating due to a singular matrix.
!split
===== Mathematical Properties =====
@@ -1753,54 +1803,6 @@ We will come back to more interpreations after we have gone through some of the
For more discussions of Ridge and Lasso regression, "Wessel van Wieringen's":"https://arxiv.org/abs/1509.09169" article is highly recommended.
Similarly, "Mehta et al's article":"https://arxiv.org/abs/1803.08823" is also recommended.
!split
===== Codes for the SVD =====
!bc pycod
import numpy as np
# SVD inversion
def SVDinv(A):
''' Takes as input a numpy matrix A and returns inv(A) based on singular value decomposition (SVD).
SVD is numerically more stable than the inversion algorithms provided by
numpy and scipy.linalg at the cost of being slower.
'''
U, s, VT = np.linalg.svd(A)
# print('test U')
# print( (np.transpose(U) @ U - U @np.transpose(U)))
# print('test VT')
# print( (np.transpose(VT) @ VT - VT @np.transpose(VT)))
print(U)
print(s)
print(VT)
D = np.zeros((len(U),len(VT)))
for i in range(0,len(VT)):
D[i,i]=s[i]
UT = np.transpose(U); V = np.transpose(VT); invD = np.linalg.inv(D)
return np.matmul(V,np.matmul(invD,UT))
X = np.array([ [1.0, -1.0, 2.0], [1.0, 0.0, 1.0], [1.0, 2.0, -1.0], [1.0, 1.0, 0.0] ])
print(X)
A = np.transpose(X) @ X
print(A)
# Brute force inversion of super-collinear matrix
#B = np.linalg.inv(A)
#print(B)
C = SVDinv(A)
print(C)
!ec
The matrix $\bm{X}$ has columns that are linearly dependent. The first
column is the row-wise sum of the other two columns. The rank of a
matrix (the column rank) is the dimension of space spanned by the
column vectors. The rank of the matrix is the number of linearly
independent columns, in this case just $2$. We see this from the
singular values when running the above code. Running the standard
inversion algorithm for matrix inversion with $\bm{X}^T\bm{X}$ results
in the program terminating due to a singular matrix.
!split
===== A better understanding of regularization =====