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

test

Browse Source
This commit is contained in:
Maciej Czajka 2022-05-18 11:54:15 +02:00
parent 85ba70d8b3
commit bddd173633
1 changed files with 124 additions and 41 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

165
jupyter.ipynb
View File

@ -2,7 +2,7 @@
"cells": [ "cells": [
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 1, "execution_count": 163,
"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_57719/1134982733.py:12: 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_61525/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"
] ]
}, },
@ -41,6 +41,9 @@
"import scipy.linalg as la\n", "import scipy.linalg as la\n",
"from PIL import Image\n", "from PIL import Image\n",
"from ipywidgets import interact\n", "from ipywidgets import interact\n",
"from numpy.linalg import eig\n",
"from math import isclose\n",
"import pprint\n",
"\n", "\n",
"# zmień szerokość komórki\n", "# zmień szerokość komórki\n",
"from IPython.core.display import display, HTML\n", "from IPython.core.display import display, HTML\n",
@ -49,7 +52,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 2, "execution_count": 164,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -115,7 +118,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 3, "execution_count": 284,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -124,7 +127,7 @@
"outputs": [], "outputs": [],
"source": [ "source": [
"gray_images = {\n", "gray_images = {\n",
" \"cat\":rgb2gray(img_as_float(data.chelsea())),\n", " \"cat\":rgb2gray(img_as_float(data.chelsea()))/255.0,\n",
" \"astro\":rgb2gray(img_as_float(data.astronaut())),\n", " \"astro\":rgb2gray(img_as_float(data.astronaut())),\n",
" \"camera\":data.camera(),\n", " \"camera\":data.camera(),\n",
" \"coin\": data.coins(),\n", " \"coin\": data.coins(),\n",
@ -158,51 +161,113 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 4, "execution_count": 304,
"metadata": {}, "metadata": {},
"outputs": [], "outputs": [
{
"ename": "IndexError",
"evalue": "index 3 is out of bounds for axis 0 with size 3",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mIndexError\u001b[0m Traceback (most recent call last)",
"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",
"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",
"\u001b[0;31mIndexError\u001b[0m: index 3 is out of bounds for axis 0 with size 3"
]
}
],
"source": [ "source": [
"def calculate(A):\n", "def calculate(A):\n",
" m = A.shape[0]\n", " n = A.shape[0]\n",
" n = 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", " # 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",
" eigenvalues, eigenvectors = la.eigh(helper)\n", "# print(f'helper ---------- {helper}')\n",
" eigenvalues, eigenvectors = eig(helper)\n",
"# print(eigenvalues)\n",
"# print(eigenvectors)\n",
" # descending sort of all the eigenvectors according to their eigenvalues\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",
" eigenvectors = np.real(eigenvectors)\n",
" eigenvalues = eigenvalues[index]\n", " eigenvalues = eigenvalues[index]\n",
" 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",
" eigenvalues2, eigenvectors2 = eig(helper2)\n",
" # descending sort of all the eigenvectors according to their eigenvalues\n",
" index2 = eigenvalues2.argsort()[::-1]\n",
"# print(f'e2 ----------------- {np.real(eigenvalues2)}')\n",
" eigenvalues2 = np.real(eigenvalues2)\n",
" eigenvectors2 = np.real(eigenvectors2)\n",
" eigenvalues2 = eigenvalues2[index2]\n",
" eigenvectors2 = eigenvectors2[:, index2]\n",
"# print(eigenvalues2)\n",
"# print(eigenvectors2)\n",
" V = np.transpose(eigenvectors2)\n",
" \n",
" # S is a diagonal matrix that keeps square root of eigenvalues\n", " # S is a diagonal matrix that keeps square root of eigenvalues\n",
" j = 0\n", " j = 0\n",
" for i in eigenvalues:\n", " for i in eigenvalues2:\n",
" if j == S.size:\n", " if j == S.size:\n",
" break\n", " break\n",
" elif i >= 0:\n", " elif i >= 0:\n",
" S[j] = np.sqrt(i)\n", " S[j] = np.sqrt(i)\n",
" j += 1\n", " j += 1\n",
" # same finding process for transpose(A)*A\n",
" helper = np.matmul(np.transpose(A), A)\n",
" eigenvalues, eigenvectors = la.eigh(helper)\n",
" # descending sort of all the eigenvectors according to their eigenvalues\n",
" index = eigenvalues.argsort()[::-1]\n",
" eigenvalues = eigenvalues[index]\n",
" eigenvectors = eigenvectors[:, index]\n",
" V = np.transpose(eigenvectors)\n",
"\n", "\n",
" # sorting S in descending order\n", " # sorting S in descending order\n",
" S[::-1].sort()\n", " S[::-1].sort()\n",
" # print_to_file(S)\n", " # print_to_file(S)\n",
"\n", "\n",
" return U, S, V" " return U, S, V\n",
"\n",
"def check_sign_V(V, builtin):\n",
" for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n",
" if builtin[j][i] < 0.0 and V[j][i] > 0.0:\n",
" V[j][i] *= -1\n",
" elif builtin[j][i] > 0.0 and V[j][i] < 0.0:\n",
" V[j][i] *= -1\n",
"\n",
" return V\n",
"\n",
"\n",
"def check_sign_U(U, builtin):\n",
" for i in range(0,builtin.shape[0]):\n",
" for j in range(0,builtin.shape[1]):\n",
" if builtin[j][i] < 0.0 and U[j][i] > 0.0:\n",
" U[j][i] *= -1\n",
" elif builtin[j][i] > 0.0 and U[j][i] < 0.0:\n",
" U[j][i] *= -1\n",
"\n",
" return U\n",
"\n",
"\n",
"a = np.array([[1,0, 2], \n",
" [0, 0, 0],\n",
" [0, 0, 0],\n",
" [0, 4, 0]])\n",
"U,s,V = calculate(np.transpose(a))\n",
"U1,s1,V1 = svd(a)\n",
"U = check_sign_U(U, U1)\n",
"V = check_sign_V(V, V1)\n",
"print(V.shape[0])\n",
"print('U: ', U, '\\n', 's: ',s, '\\n','V: ',V, '\\n')\n",
"print('--------------------------------------')\n",
"print('U1: ',U1, '\\n','s1: ',s1, '\\n', 'V1: ',V1, '\\n')\n",
"k = 4\n",
"reconst_matrix4 = np.dot(U[:,:k],np.dot(np.diag(s[:k]),V[:k,:]))\n",
"print('reconst -------------- ', reconst_matrix4)"
] ]
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 5, "execution_count": 293,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -214,11 +279,27 @@
" \"\"\"\n", " \"\"\"\n",
" 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",
"# U,s,V = svd(image,full_matrices=False)\n", " image = np.real(image)\n",
" U,s,V = calculate(image)\n", " print(f'image ---------- {image}')\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_matrix,s" " return reconst_matrix2,s"
] ]
}, },
{ {
@ -235,7 +316,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 6, "execution_count": 294,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -274,7 +355,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 7, "execution_count": 295,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -304,7 +385,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 21, "execution_count": 296,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -324,8 +405,8 @@
" original_shape = image.shape\n", " original_shape = image.shape\n",
" 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", "# 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",
@ -338,7 +419,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 22, "execution_count": 297,
"metadata": { "metadata": {
"pycharm": { "pycharm": {
"name": "#%%\n" "name": "#%%\n"
@ -348,12 +429,12 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "7c4a40e6f857407586d21ab113128a97", "model_id": "092923713cd74de2b08c173639aae9d4",
"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', 'coffee', 'rocket', 'koala', '…" "interactive(children=(Dropdown(description='img_name', options=('cat', 'astro', 'camera', 'coin', 'clock', 'te…"
] ]
}, },
"metadata": {}, "metadata": {},
@ -365,7 +446,7 @@
"<function __main__.print_matrices(img_name, k)>" "<function __main__.print_matrices(img_name, k)>"
] ]
}, },
"execution_count": 22, "execution_count": 297,
"metadata": {}, "metadata": {},
"output_type": "execute_result" "output_type": "execute_result"
} }
@ -376,17 +457,18 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 10, "execution_count": 298,
"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": "3533923ce918489eab6fc7de25234078", "model_id": "95f52e9c633540678b6285ce25adaff0",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -551,13 +633,14 @@
"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": "0af20b426d1240cca1e08e40ede10e98", "model_id": "9976c068c2cc40788cdce8bf0365c66a",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -595,7 +678,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "86f500d175554478ac9805075a809d28", "model_id": "24661e5fa1644a349470561f5d9b1324",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -703,7 +786,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "e97ca096507f42b49f98dbb87ac5a2ec", "model_id": "3026b3e43b6c4921abd17c4cf37e9954",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },
@ -731,7 +814,7 @@
{ {
"data": { "data": {
"application/vnd.jupyter.widget-view+json": { "application/vnd.jupyter.widget-view+json": {
"model_id": "71fea739144848d69afb1662172190ad", "model_id": "a93634f89fdf459189067b59638f3021",
"version_major": 2, "version_major": 2,
"version_minor": 0 "version_minor": 0
}, },