s444018/svd_mpsic
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

final version

Browse Source
This commit is contained in:
Maciej Czajka 2022-05-18 12:44:39 +02:00
parent bddd173633
commit 1bde756004
2 changed files with 89 additions and 87 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

152
jupyter.ipynb
View File

@ -2,7 +2,7 @@
"cells": [ "cells": [
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 163, "execution_count": 1,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"is_executing": true, "is_executing": true,
@ -14,7 +14,7 @@
"name": "stderr", "name": "stderr",
"output_type": "stream", "output_type": "stream",
"text": [ "text": [
"/var/folders/tq/jq5nwbnj7v10tls99x99qbh40000gn/T/ipykernel_61525/3026179557.py:15: DeprecationWarning: Importing display from IPython.core.display is deprecated since IPython 7.14, please import from IPython display\n", "/var/folders/tq/jq5nwbnj7v10tls99x99qbh40000gn/T/ipykernel_65173/3026179557.py:15: DeprecationWarning: Importing display from IPython.core.display is deprecated since IPython 7.14, please import from IPython display\n",
" from IPython.core.display import display, HTML\n" " from IPython.core.display import display, HTML\n"
] ]
}, },
@ -52,7 +52,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 164, "execution_count": 2,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -118,7 +118,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 284, "execution_count": 3,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -127,8 +127,9 @@
"outputs": [], "outputs": [],
"source": [ "source": [
"gray_images = {\n", "gray_images = {\n",
" \"cat\":rgb2gray(img_as_float(data.chelsea()))/255.0,\n", " \"cat\":rgb2gray(img_as_float(data.chelsea())),\n",
" \"astro\":rgb2gray(img_as_float(data.astronaut())),\n", " \"astro\":rgb2gray(img_as_float(data.astronaut())),\n",
" \"line\": rgb2gray(img_as_float(Image.open('line.jpeg'))),\n",
" \"camera\":data.camera(),\n", " \"camera\":data.camera(),\n",
" \"coin\": data.coins(),\n", " \"coin\": data.coins(),\n",
" \"clock\":data.clock(),\n", " \"clock\":data.clock(),\n",
@ -161,19 +162,41 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 304, "execution_count": 4,
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
{ {
"ename": "IndexError", "name": "stdout",
"evalue": "index 3 is out of bounds for axis 0 with size 3", "output_type": "stream",
"output_type": "error", "text": [
"traceback": [ "5\n",
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "U: [[ 0. 0. 1. 0.]\n",
"\u001b[0;31mIndexError\u001b[0m Traceback (most recent call last)", " [ 0. 1. 0. 0.]\n",
"Input \u001b[0;32mIn [304]\u001b[0m, in \u001b[0;36m<cell line: 77>\u001b[0;34m()\u001b[0m\n\u001b[1;32m 75\u001b[0m U,s,V \u001b[38;5;241m=\u001b[39m calculate(np\u001b[38;5;241m.\u001b[39mtranspose(a))\n\u001b[1;32m 76\u001b[0m U1,s1,V1 \u001b[38;5;241m=\u001b[39m svd(a)\n\u001b[0;32m---> 77\u001b[0m U \u001b[38;5;241m=\u001b[39m \u001b[43mcheck_sign_U\u001b[49m\u001b[43m(\u001b[49m\u001b[43mU\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mU1\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 78\u001b[0m V \u001b[38;5;241m=\u001b[39m check_sign_V(V, V1)\n\u001b[1;32m 79\u001b[0m \u001b[38;5;28mprint\u001b[39m(V\u001b[38;5;241m.\u001b[39mshape[\u001b[38;5;241m0\u001b[39m])\n", " [ 0. 0. 0. -1.]\n",
"Input \u001b[0;32mIn [304]\u001b[0m, in \u001b[0;36mcheck_sign_U\u001b[0;34m(U, builtin)\u001b[0m\n\u001b[1;32m 61\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m i \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28mrange\u001b[39m(\u001b[38;5;241m0\u001b[39m,builtin\u001b[38;5;241m.\u001b[39mshape[\u001b[38;5;241m0\u001b[39m]):\n\u001b[1;32m 62\u001b[0m \u001b[38;5;28;01mfor\u001b[39;00m j \u001b[38;5;129;01min\u001b[39;00m \u001b[38;5;28mrange\u001b[39m(\u001b[38;5;241m0\u001b[39m,builtin\u001b[38;5;241m.\u001b[39mshape[\u001b[38;5;241m1\u001b[39m]):\n\u001b[0;32m---> 63\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m builtin[j][i] \u001b[38;5;241m<\u001b[39m \u001b[38;5;241m0.0\u001b[39m \u001b[38;5;129;01mand\u001b[39;00m \u001b[43mU\u001b[49m\u001b[43m[\u001b[49m\u001b[43mj\u001b[49m\u001b[43m]\u001b[49m[i] \u001b[38;5;241m>\u001b[39m \u001b[38;5;241m0.0\u001b[39m:\n\u001b[1;32m 64\u001b[0m U[j][i] \u001b[38;5;241m*\u001b[39m\u001b[38;5;241m=\u001b[39m \u001b[38;5;241m-\u001b[39m\u001b[38;5;241m1\u001b[39m\n\u001b[1;32m 65\u001b[0m \u001b[38;5;28;01melif\u001b[39;00m builtin[j][i] \u001b[38;5;241m>\u001b[39m \u001b[38;5;241m0.0\u001b[39m \u001b[38;5;129;01mand\u001b[39;00m U[j][i] \u001b[38;5;241m<\u001b[39m \u001b[38;5;241m0.0\u001b[39m:\n", " [ 1. 0. 0. 0.]] \n",
"\u001b[0;31mIndexError\u001b[0m: index 3 is out of bounds for axis 0 with size 3" " s: [4. 3. 2.23606798 0. ] \n",
" V: [[ 0. 1. 0. 0. 0. ]\n",
" [ 0. 0. 1. 0. 0. ]\n",
" [ 0.4472136 0. 0. 0. 0.89442719]\n",
" [ 0. 0. 0. 1. 0. ]\n",
" [-0.89442719 0. 0. 0. 0.4472136 ]] \n",
"\n",
"--------------------------------------\n",
"U1: [[ 0. 0. 1. 0.]\n",
" [ 0. 1. 0. 0.]\n",
" [ 0. 0. 0. -1.]\n",
" [ 1. 0. 0. 0.]] \n",
" s1: [4. 3. 2.23606798 0. ] \n",
" V1: [[-0. 1. 0. 0. 0. ]\n",
" [-0. 0. 1. 0. 0. ]\n",
" [ 0.4472136 0. 0. 0. 0.89442719]\n",
" [ 0. 0. 0. 1. 0. ]\n",
" [-0.89442719 0. 0. 0. 0.4472136 ]] \n",
"\n",
"reconst -------------- [[1. 0. 0. 0. 2.]\n",
" [0. 0. 3. 0. 0.]\n",
" [0. 0. 0. 0. 0.]\n",
" [0. 4. 0. 0. 0.]]\n"
] ]
} }
], ],
@ -183,13 +206,9 @@
" m = A.shape[1]\n", " m = A.shape[1]\n",
" S = np.zeros(n)\n", " S = np.zeros(n)\n",
"\n", "\n",
" # finding eigenvectors with biggest eigenvalues of A*transpose(A)\n",
" helper = np.matmul(A, np.transpose(A))\n", " helper = np.matmul(A, np.transpose(A))\n",
"# print(f'helper ---------- {helper}')\n",
" eigenvalues, eigenvectors = eig(helper)\n", " eigenvalues, eigenvectors = eig(helper)\n",
"# print(eigenvalues)\n", " \n",
"# print(eigenvectors)\n",
" # descending sort of all the eigenvectors according to their eigenvalues\n",
" index = eigenvalues.argsort()[::-1]\n", " index = eigenvalues.argsort()[::-1]\n",
" eigenvalues = np.real(eigenvalues)\n", " eigenvalues = np.real(eigenvalues)\n",
" eigenvectors = np.real(eigenvectors)\n", " eigenvectors = np.real(eigenvectors)\n",
@ -197,21 +216,15 @@
" eigenvectors = eigenvectors[:, index]\n", " eigenvectors = eigenvectors[:, index]\n",
" U = eigenvectors\n", " U = eigenvectors\n",
"\n", "\n",
" # same finding process for transpose(A)*A\n",
" helper2 = np.matmul(np.transpose(A), A)\n", " helper2 = np.matmul(np.transpose(A), A)\n",
" eigenvalues2, eigenvectors2 = eig(helper2)\n", " eigenvalues2, eigenvectors2 = eig(helper2)\n",
" # descending sort of all the eigenvectors according to their eigenvalues\n",
" index2 = eigenvalues2.argsort()[::-1]\n", " index2 = eigenvalues2.argsort()[::-1]\n",
"# print(f'e2 ----------------- {np.real(eigenvalues2)}')\n",
" eigenvalues2 = np.real(eigenvalues2)\n", " eigenvalues2 = np.real(eigenvalues2)\n",
" eigenvectors2 = np.real(eigenvectors2)\n", " eigenvectors2 = np.real(eigenvectors2)\n",
" eigenvalues2 = eigenvalues2[index2]\n", " eigenvalues2 = eigenvalues2[index2]\n",
" eigenvectors2 = eigenvectors2[:, index2]\n", " eigenvectors2 = eigenvectors2[:, index2]\n",
"# print(eigenvalues2)\n",
"# print(eigenvectors2)\n",
" V = np.transpose(eigenvectors2)\n", " V = np.transpose(eigenvectors2)\n",
" \n", " \n",
" # S is a diagonal matrix that keeps square root of eigenvalues\n",
" j = 0\n", " j = 0\n",
" for i in eigenvalues2:\n", " for i in eigenvalues2:\n",
" if j == S.size:\n", " if j == S.size:\n",
@ -220,39 +233,53 @@
" S[j] = np.sqrt(i)\n", " S[j] = np.sqrt(i)\n",
" j += 1\n", " j += 1\n",
"\n", "\n",
" # sorting S in descending order\n",
" S[::-1].sort()\n", " S[::-1].sort()\n",
" # print_to_file(S)\n",
"\n", "\n",
" return U, S, V\n", " return U, S, V\n",
"\n", "\n",
"def check_sign_V(V, builtin):\n", "def check_sign_V(V, builtin):\n",
" if builtin.shape[0] < builtin.shape[1]:\n",
" for i in range(0,builtin.shape[0]):\n", " for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n", " for j in range(0,builtin.shape[1]):\n",
" if builtin[j][i] < 0.0 and V[j][i] > 0.0:\n", " if builtin[j][i] < 0.0 and V[j][i] > 0.0:\n",
" V[j][i] *= -1\n", " V[j][i] *= -1\n",
" elif builtin[j][i] > 0.0 and V[j][i] < 0.0:\n", " elif builtin[j][i] > 0.0 and V[j][i] < 0.0:\n",
" V[j][i] *= -1\n", " V[j][i] *= -1\n",
" else:\n",
" for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n",
" if builtin[i][j] < 0.0 and V[i][j] > 0.0:\n",
" V[i][j] *= -1\n",
" elif builtin[i][j] > 0.0 and V[i][j] < 0.0:\n",
" V[i][j] *= -1\n",
"\n", "\n",
" return V\n", " return V\n",
"\n", "\n",
"\n", "\n",
"def check_sign_U(U, builtin):\n", "def check_sign_U(U, builtin):\n",
" if builtin.shape[0] < builtin.shape[1]: \n",
" for i in range(0,builtin.shape[0]):\n", " for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n", " for j in range(0,builtin.shape[1]):\n",
" if builtin[j][i] < 0.0 and U[j][i] > 0.0:\n", " if builtin[j][i] < 0.0 and U[j][i] > 0.0:\n",
" U[j][i] *= -1\n", " U[j][i] *= -1\n",
" elif builtin[j][i] > 0.0 and U[j][i] < 0.0:\n", " elif builtin[j][i] > 0.0 and U[j][i] < 0.0:\n",
" U[j][i] *= -1\n", " U[j][i] *= -1\n",
" else:\n",
" for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n",
" if builtin[i][j] < 0.0 and U[i][j] > 0.0:\n",
" U[i][j] *= -1\n",
" elif builtin[j][j] > 0.0 and U[i][j] < 0.0:\n",
" U[i][j] *= -1\n",
"\n", "\n",
" return U\n", " return U\n",
"\n", "\n",
"\n", "\n",
"a = np.array([[1,0, 2], \n", "a = np.array([[1,0, 0, 0, 2], \n",
" [0, 0, 0],\n", " [0, 0, 3, 0, 0],\n",
" [0, 0, 0],\n", " [0, 0, 0, 0, 0],\n",
" [0, 4, 0]])\n", " [0, 4, 0, 0, 0]])\n",
"U,s,V = calculate(np.transpose(a))\n", "U,s,V = calculate(a)\n",
"U1,s1,V1 = svd(a)\n", "U1,s1,V1 = svd(a)\n",
"U = check_sign_U(U, U1)\n", "U = check_sign_U(U, U1)\n",
"V = check_sign_V(V, V1)\n", "V = check_sign_V(V, V1)\n",
@ -267,7 +294,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 293, "execution_count": 5,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -280,26 +307,10 @@
" Wykonaj dekompozycję SVD, a następnie okrojoną rekonstrukcję (przy użyciu k wartości/wektorów osobliwych)\n", " Wykonaj dekompozycję SVD, a następnie okrojoną rekonstrukcję (przy użyciu k wartości/wektorów osobliwych)\n",
" \"\"\"\n", " \"\"\"\n",
" image = np.real(image)\n", " image = np.real(image)\n",
" print(f'image ---------- {image}')\n",
" U,s,V = svd(image,full_matrices=False)\n", " U,s,V = svd(image,full_matrices=False)\n",
" U1,s1,V1 = calculate(image)\n",
" U1 = np.real(U1)\n",
" s1 = np.real(s1)\n",
" V1 = np.real(V1)\n",
"# U1 = check_sign_U(U1, U)\n",
"# V1 = check_sign_V(V1, V)\n",
" print(U.shape, U1.shape)\n",
" print(V.shape, V1.shape)\n",
" print(s.shape, s1.shape)\n",
" reconst_matrix = np.dot(U[:,:k],np.dot(np.diag(s[:k]),V[:k,:]))\n", " reconst_matrix = np.dot(U[:,:k],np.dot(np.diag(s[:k]),V[:k,:]))\n",
" S = np.zeros(n)\n",
" for i in range(0,k):\n",
" suma[:,:i] += np.dot(np.dot(s1[i],U1[:,:i]), V1[:,:i])\n",
" print(suma)\n",
" reconst_matrix2 = np.dot(np.dot(np.diag(s1[:k]),U1[:,:k]),V1[:,:k])\n",
" print(reconst_matrix, '============================================', reconst_matrix2)\n",
"\n", "\n",
" return reconst_matrix2,s" " return reconst_matrix,s"
] ]
}, },
{ {
@ -316,7 +327,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 294, "execution_count": 6,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -355,7 +366,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 295, "execution_count": 7,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -385,7 +396,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 296, "execution_count": 8,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -406,7 +417,6 @@
" print(f\"Input image dimensions. Width:{original_shape[1]} Height:{original_shape[0]}\")\n", " print(f\"Input image dimensions. Width:{original_shape[1]} Height:{original_shape[0]}\")\n",
"\n", "\n",
" U,s,V = svd(image,full_matrices=False)\n", " U,s,V = svd(image,full_matrices=False)\n",
"# U,s,V = calculate(image)\n",
" print(f\"Shape of U matrix: {U[:,:k].shape}\")\n", " print(f\"Shape of U matrix: {U[:,:k].shape}\")\n",
" print(f\"U MATRIX: {U[:,:k]}\")\n", " print(f\"U MATRIX: {U[:,:k]}\")\n",
" print('*' * 100)\n", " print('*' * 100)\n",
@ -419,7 +429,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 297, "execution_count": 9,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -429,12 +439,12 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "092923713cd74de2b08c173639aae9d4", "model_id": "702fd8dcaa0c495d8944c5b0dd82892b",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
"text/plain": [ "text/plain": [
"interactive(children=(Dropdown(description='img_name', options=('cat', 'astro', 'camera', 'coin', 'clock', 'te…" "interactive(children=(Dropdown(description='img_name', options=('cat', 'astro', 'line', 'camera', 'coin', 'clo…"
] ]
}, },
"metadata": {}, "metadata": {},
@ -446,7 +456,7 @@
"<function __main__.print_matrices(img_name, k)>" "<function __main__.print_matrices(img_name, k)>"
] ]
}, },
"execution_count": 297, "execution_count": 9,
"metadata": {}, "metadata": {},
"output_type": "execute_result" "output_type": "execute_result"
} }
@ -457,7 +467,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 298, "execution_count": 10,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -468,12 +478,12 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "95f52e9c633540678b6285ce25adaff0", "model_id": "405c92f2998843498e9cd20675ba2b4e",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
"text/plain": [ "text/plain": [
"interactive(children=(Dropdown(description='img_name', options=('cat', 'astro', 'camera', 'coin', 'clock', 'te…" "interactive(children=(Dropdown(description='img_name', options=('cat', 'astro', 'line', 'camera', 'coin', 'clo…"
] ]
}, },
"metadata": {}, "metadata": {},
@ -633,14 +643,13 @@
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
}, }
"scrolled": false
}, },
"outputs": [ "outputs": [
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "9976c068c2cc40788cdce8bf0365c66a", "model_id": "6d2c561f500d42be903e30933d3fe7d6",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -678,7 +687,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "24661e5fa1644a349470561f5d9b1324", "model_id": "f4b99140de5d495580f6abaa20a04b4c",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -786,7 +795,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "3026b3e43b6c4921abd17c4cf37e9954", "model_id": "6363d12d2fbe49518bf5735f6472a497",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -814,7 +823,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "a93634f89fdf459189067b59638f3021", "model_id": "edd32a03e6184b15ae3283c8e24f2d3e",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -829,13 +838,6 @@
"source": [ "source": [
"interact(compress_show_color_images_layer,img_name=list_widget,k=int_slider_widget);" "interact(compress_show_color_images_layer,img_name=list_widget,k=int_slider_widget);"
] ]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
} }
], ],
"metadata": { "metadata": {

BIN
line.jpeg Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 8.1 KiB