** Algorithmes d'optimisation -- L3 MINT et doubles licences 2017/2018 -- Université Paris-Sud **
$\newcommand{\Rsp}{\mathbb{R}} \newcommand{\nr}[1]{\|#1\|} \newcommand{\abs}[1]{|#1|} \newcommand{\eps}{\varepsilon} \newcommand{\sca}[2]{\langle#1|#2\rangle} \newcommand{\D}{\mathrm{D}} \newcommand{\hdots}{\dots} \newcommand{\cond}{\mathrm{cond}}$
L'explication de la méthode est à lire dans la feuille de TD. On commence par construire les données $(x_i)_{i\in I \cup J} \subseteq \Rsp^2$:
On cherche alors à construire, par optimisation de la fonction $F$, un vecteur $w\in\Rsp^2$ tel que $\sigma(\sca{w}{x_i}) \simeq 0$ si $i\in I$ et $\sigma(\sca{w}{x_i})\simeq 1$ si $i\in J$.
On définit également des fonctions tirées des TP précédents:
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
# la commande suivante agrandit les figures
plt.rcParams['figure.figsize'] = [9.,6.]
# initialiser les variables K,Q, x0, fmin
n = 30
X = np.vstack((np.hstack((-0.5+.2*np.random.randn(n,1), -0.5 + .2*np.random.randn(n,1))),
np.hstack((0.5+.2*np.random.randn(n,1), 0.5 + .2*np.random.randn(n,1)))))
I = range(0,n)
J = range(n,2*n)
plt.plot(X[I,0],X[I,1],'.r')
plt.plot(X[J,0],X[J,1],'.g')
plt.axis('equal');
# on recupere ici des fonctions issues des TP précédents. Seule la fonction check_hessian est nouvelle.
def backtrack(f,x,d,m,alpha=0.3,beta=0.5):
t = 1
while f(x+t*d) > f(x) + alpha*t*m:
t = beta*t
return t
def gradient_backtracking(f,g,x0,err=1e-6,maxiter=500):
x = x0.copy()
fiter = []
giter = []
k = 0 # nombre d'itérations
while(True):
k = k+1
if k > maxiter: # maximum de 10^6 itérations
print('erreur: nombre maximum d\'itérations atteint')
break
d = -g(x)
fiter.append(f(x))
giter.append(np.linalg.norm(d))
if np.linalg.norm(d) <= err:
break
t = backtrack(f,x,d,-np.linalg.norm(d)**2)
if k%10==0: # on affiche des informations toute les 20 itérations
print('iteration %d: f=%f, |g|=%f, step=%f' % (k, f(x), np.linalg.norm(d),t))
x = x + t*d
return x,np.array(fiter),np.array(giter)
def check_gradient(f,g,x0):
N = len(x0)
gg = np.zeros(N)
for i in range(N):
eps = 1e-4
e = np.zeros(N)
e[i] = eps
gg[i] = (f(x0+e) - f(x0-e))/(2*eps)
print('erreur numerique dans le calcul du gradient: %g (doit etre petit)' % np.linalg.norm(g(x0)-gg))
# cette fonction vérifie que h (=fonction évaluant le hessian) est bien la dérivée
# de g (= fonction évaluant le gradient), en un point x0
def check_hessian(g,h,x0):
N = len(x0)
H = np.zeros((N,N))
for i in range(N):
eps = 1e-5
e = np.zeros(N)
e[i] = eps
H[i,:] = (g(x0+e) - g(x0-e))/(2*eps)
print('erreur numerique dans le calcul de la hessienne: %g (doit etre petit)' % np.sum((H-h(x0)))**2)
QI.1: Écrire une fonction $F(w,X,I,J,\gamma)$ et $gF(w,X,I,J,\gamma)$ calculant la valeur et le gradient de la fonctionnelle
$$ F(w) = - \left(\sum_{i\in I} \log(1-\sigma(\sca{w}{x_i})) + \sum_{i\in J} \log(\sigma(\sca{w}{x_i})\right) + \frac{\gamma}{2} \nr{w}^2. $$Tester que $gF$ correspond bien au gradient de $F$ en utilisant la fonction check_gradient.
# on définit une fonction sigma(t) correspondant à la sigmoïde
def sigma(t):
return 1/(1+np.exp(-t))
def F(w,X,I,J,gamma):
# on calcule F en utilisant la formule (on pourra utiliser np.sum)
h = sigma(np.dot(X,w.T))
return -(np.sum(np.log(1.-h[I])) + np.sum(np.log(h[J]))) + .5*gamma*np.dot(w,w)
def gF(w,X,I,J,gamma):
d = len(w)
g = np.zeros(d)
# on calcule le gradient en appliquant la formule obtenue en TD
h = sigma(np.dot(X,w.T))
for i in I: g += h[i]*X[i,:]
for j in J: g += (h[j]-1)*X[j,:]
g += gamma*w
return g
# on vérifie le calcul du gradient en utilisant la fonction check_gradient
# rappel: on peut utiliser f = lambda w: F(w,X,I,J,gamma) pour définir
# f(w) = F(w,X,I,J,gamma)
gamma = .1
check_gradient(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma), np.random.rand(2))
erreur numerique dans le calcul du gradient: 2.4316e-09 (doit etre petit)
QI.2: Calculer $w=\arg\min_{\Rsp^2} F$ en utilisant la fonction gradient_backtracking. Tracer sur un graphe l'évolution de la norme du gradient au cours des itérations de la descente de gradient, soit $(k,\nr{\nabla F(w^{(k)})})$, où l'axe $y$ sera en échelle logarithmique.
w0 = np.zeros(2)
gamma = .01
w,fiter,giter = gradient_backtracking(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma), w0, err=1e-7, maxiter=1000)
plt.semilogy(giter)
iteration 10: f=0.869162, |g|=0.138801, step=1.000000 iteration 20: f=0.783502, |g|=0.054698, step=1.000000 iteration 30: f=0.766904, |g|=0.027903, step=1.000000 iteration 40: f=0.762141, |g|=0.016022, step=1.000000 iteration 50: f=0.760460, |g|=0.009998, step=1.000000 iteration 60: f=0.759767, |g|=0.006669, step=1.000000 iteration 70: f=0.759445, |g|=0.004682, step=1.000000 iteration 80: f=0.759281, |g|=0.003406, step=1.000000 iteration 90: f=0.759192, |g|=0.002534, step=1.000000 iteration 100: f=0.759143, |g|=0.001910, step=1.000000 iteration 110: f=0.759114, |g|=0.001450, step=1.000000 iteration 120: f=0.759098, |g|=0.001106, step=1.000000 iteration 130: f=0.759088, |g|=0.000846, step=1.000000 iteration 140: f=0.759082, |g|=0.000648, step=1.000000 iteration 150: f=0.759079, |g|=0.000496, step=1.000000 iteration 160: f=0.759077, |g|=0.000380, step=1.000000 iteration 170: f=0.759076, |g|=0.000292, step=1.000000 iteration 180: f=0.759075, |g|=0.000224, step=1.000000 iteration 190: f=0.759075, |g|=0.000171, step=1.000000 iteration 200: f=0.759075, |g|=0.000132, step=1.000000 iteration 210: f=0.759075, |g|=0.000101, step=1.000000 iteration 220: f=0.759075, |g|=0.000077, step=1.000000 iteration 230: f=0.759075, |g|=0.000059, step=1.000000 iteration 240: f=0.759075, |g|=0.000046, step=1.000000 iteration 250: f=0.759075, |g|=0.000035, step=1.000000 iteration 260: f=0.759075, |g|=0.000027, step=1.000000 iteration 270: f=0.759075, |g|=0.000021, step=1.000000 iteration 280: f=0.759074, |g|=0.000016, step=1.000000 iteration 290: f=0.759074, |g|=0.000012, step=1.000000 iteration 300: f=0.759074, |g|=0.000009, step=1.000000 iteration 310: f=0.759074, |g|=0.000007, step=1.000000 iteration 320: f=0.759074, |g|=0.000005, step=1.000000 iteration 330: f=0.759074, |g|=0.000004, step=1.000000 iteration 340: f=0.759074, |g|=0.000003, step=1.000000 iteration 350: f=0.759074, |g|=0.000002, step=1.000000 iteration 360: f=0.759074, |g|=0.000002, step=1.000000 iteration 370: f=0.759074, |g|=0.000001, step=1.000000 iteration 380: f=0.759074, |g|=0.000001, step=1.000000 iteration 390: f=0.759074, |g|=0.000001, step=1.000000 iteration 400: f=0.759074, |g|=0.000001, step=1.000000 iteration 410: f=0.759074, |g|=0.000001, step=1.000000 iteration 420: f=0.759074, |g|=0.000000, step=1.000000 iteration 430: f=0.759074, |g|=0.000000, step=1.000000 iteration 440: f=0.759074, |g|=0.000000, step=1.000000 iteration 450: f=0.759074, |g|=0.000000, step=1.000000 iteration 460: f=0.759074, |g|=0.000000, step=1.000000 iteration 470: f=0.759074, |g|=0.000000, step=1.000000
[<matplotlib.lines.Line2D at 0x7f80fd446e50>]
QI.3: Soit $w$ le minimiseur approché de $F$ calculé dans la question précédente. Tracer sur un graphe:
Faire varier le paramètre de régularisation $\gamma \in \{0.01,0.1,1,10\}$ et observer:
# on définit une fonction permettant d'afficher les lignes de niveaux de u: R^2 -> R
# NB: u peut être définie en utilisant lambda, par exemple:
# u = lambda (x,y): ... où x,y sont les coordonnées du point
def lignes_niveau(u):
Xcontour,Ycontour = np.meshgrid(np.linspace(-1., 1., 100),
np.linspace(-1., 1., 100))
Zcontour = u(Xcontour,Ycontour)
plt.axes([-1.5, -1.5, 1.5, 1.5])
plt.axis('equal')
p = plt.contour(Xcontour, Ycontour, Zcontour, cmap='RdBu')
plt.clabel(p, inline=1, fontsize=10)
for gamma in [0.01,0.1,1,10]:
w,fiter,giter = gradient_backtracking(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma),
np.zeros(2), err=1e-7, maxiter=1000)
print('gamma = %d: niter = %d, |w|=%g' % (gamma,len(fiter),np.linalg.norm(w)))
Xcontour,Ycontour = np.meshgrid(np.linspace(-1., 1., 100),
np.linspace(-1., 1., 100))
plt.figure()
lignes_niveau(lambda x,y: sigma(x*w[0] + y*w[1]))
plt.scatter(X[I,0],X[I,1],marker='s',color='r', label='(x_j), j\in I')
plt.scatter(X[J,0],X[J,1],marker='o',color='b', label='(x_i), i\in J')
plt.axis('tight')
plt.title('donn\\\'ees $(x_i)_{i\in I\cup J}$ et solution $u: x \\mapsto \sigma(\\langle x | w\\rangle)$ / ')
plt.legend()
iteration 10: f=0.869162, |g|=0.138801, step=1.000000 iteration 20: f=0.783502, |g|=0.054698, step=1.000000 iteration 30: f=0.766904, |g|=0.027903, step=1.000000 iteration 40: f=0.762141, |g|=0.016022, step=1.000000 iteration 50: f=0.760460, |g|=0.009998, step=1.000000 iteration 60: f=0.759767, |g|=0.006669, step=1.000000 iteration 70: f=0.759445, |g|=0.004682, step=1.000000 iteration 80: f=0.759281, |g|=0.003406, step=1.000000 iteration 90: f=0.759192, |g|=0.002534, step=1.000000 iteration 100: f=0.759143, |g|=0.001910, step=1.000000 iteration 110: f=0.759114, |g|=0.001450, step=1.000000 iteration 120: f=0.759098, |g|=0.001106, step=1.000000 iteration 130: f=0.759088, |g|=0.000846, step=1.000000 iteration 140: f=0.759082, |g|=0.000648, step=1.000000 iteration 150: f=0.759079, |g|=0.000496, step=1.000000 iteration 160: f=0.759077, |g|=0.000380, step=1.000000 iteration 170: f=0.759076, |g|=0.000292, step=1.000000 iteration 180: f=0.759075, |g|=0.000224, step=1.000000 iteration 190: f=0.759075, |g|=0.000171, step=1.000000 iteration 200: f=0.759075, |g|=0.000132, step=1.000000 iteration 210: f=0.759075, |g|=0.000101, step=1.000000 iteration 220: f=0.759075, |g|=0.000077, step=1.000000 iteration 230: f=0.759075, |g|=0.000059, step=1.000000 iteration 240: f=0.759075, |g|=0.000046, step=1.000000 iteration 250: f=0.759075, |g|=0.000035, step=1.000000 iteration 260: f=0.759075, |g|=0.000027, step=1.000000 iteration 270: f=0.759075, |g|=0.000021, step=1.000000 iteration 280: f=0.759074, |g|=0.000016, step=1.000000 iteration 290: f=0.759074, |g|=0.000012, step=1.000000 iteration 300: f=0.759074, |g|=0.000009, step=1.000000 iteration 310: f=0.759074, |g|=0.000007, step=1.000000 iteration 320: f=0.759074, |g|=0.000005, step=1.000000 iteration 330: f=0.759074, |g|=0.000004, step=1.000000 iteration 340: f=0.759074, |g|=0.000003, step=1.000000 iteration 350: f=0.759074, |g|=0.000002, step=1.000000 iteration 360: f=0.759074, |g|=0.000002, step=1.000000 iteration 370: f=0.759074, |g|=0.000001, step=1.000000 iteration 380: f=0.759074, |g|=0.000001, step=1.000000 iteration 390: f=0.759074, |g|=0.000001, step=1.000000 iteration 400: f=0.759074, |g|=0.000001, step=1.000000 iteration 410: f=0.759074, |g|=0.000001, step=1.000000 iteration 420: f=0.759074, |g|=0.000000, step=1.000000 iteration 430: f=0.759074, |g|=0.000000, step=1.000000 iteration 440: f=0.759074, |g|=0.000000, step=1.000000 iteration 450: f=0.759074, |g|=0.000000, step=1.000000 iteration 460: f=0.759074, |g|=0.000000, step=1.000000 iteration 470: f=0.759074, |g|=0.000000, step=1.000000 gamma = 0: niter = 471, |w|=10.1508 iteration 10: f=3.345371, |g|=0.015735, step=1.000000 iteration 20: f=3.344646, |g|=0.002620, step=1.000000 iteration 30: f=3.344626, |g|=0.000444, step=1.000000 iteration 40: f=3.344625, |g|=0.000075, step=1.000000 iteration 50: f=3.344625, |g|=0.000013, step=1.000000 iteration 60: f=3.344625, |g|=0.000002, step=1.000000 iteration 70: f=3.344625, |g|=0.000000, step=1.000000 gamma = 0: niter = 78, |w|=6.36604 iteration 10: f=11.948578, |g|=0.000565, step=0.500000 iteration 20: f=11.948578, |g|=0.000000, step=0.500000 gamma = 1: niter = 21, |w|=3.41189 iteration 10: f=28.820387, |g|=0.000004, step=0.125000 gamma = 10: niter = 14, |w|=1.21634
QI.4: Écrire une fonction $hF$ calculant la hessienne de $F$, en utilisant la formule du TD. Vérifier son bon fonctionnement via check_hessian (voir l'introduction du TP).
def hF(w,X,I,J,K):
h = sigma(np.dot(X,w.T))
d = len(w)
N = X.shape[0]
A = K*np.eye(d)
for i in range(N):
A += h[i] * (1-h[i]) * np.tensordot(X[i,:],X[i,:],axes=0)
return A
check_hessian(lambda w: gF(w,X,I,J,gamma),
lambda w: hF(w,X,I,J,gamma), np.random.rand(2))
erreur numerique dans le calcul de la hessienne: 8.3674e-21 (doit etre petit)
Q1.5: Écrire une fonction newton_backtracking implémentant la méthode de Newton avec backtracking (ne pas hésiter à copier-coller la fonction gradient_backtracking et à changer ce qui est nécessaire) :
$$ \begin{cases} g^{(k)} = \nabla F(x^{(k)}) \\ d^{(k)} = -D^2 F(x^{(k)})^{-1} \nabla F(x^{(k)}) \\ m = \sca{d^{(k)}}{g^{(k)}} \\ t^{(k)} = \hbox{obtenu par backtracking d'Armijo, cf la fonction backtrack}\\ x^{(k+1)} = x^{(k)} + t^{(k)} d^{(k)} \end{cases} $$Utiliser cette fonction pour calculer le minimum de $F$. Comme dans la question Q1.2, tracer la norme du gradient en fonction de l'itération. Vérifier également que le nombre d'itérations dépend peu du choix du paramètre $\gamma$.
def newton_backtracking(f,g,h,x0,err=1e-6,maxiter=500):
# la boucle est tres similaire à celle de gradient_backtracking: la seule chose qui
# change est le choix de la direction de descente.
# NB: on utilise np.linalg.solve pour résoudre le système linéaire apparaissant dans la méthode
x = x0.copy()
fIter = []
gIter = []
k = 0
while(True):
k = k+1
if k > maxiter: # maximum de 10^6 itérations
print('erreur: nombre maximum d\'itérations atteint')
break
gx = -g(x)
dx = np.linalg.solve(h(x),gx)
fIter.append(f(x))
gIter.append(np.linalg.norm(dx))
if np.abs(np.linalg.norm(gx)) <= err:
break
t = backtrack(f,x,dx,np.dot(gx,dx))
print('iteration %d: f=%f, |g|=%f, step=%f' % (k, f(x), np.linalg.norm(gx),t))
x = x + t*dx
return x,np.array(fIter),np.array(gIter)
# comparer le nombre d'iterations dans la methode de Newton pour gamma dans [0.01,0.1,1,10]
for gamma in [0.01,0.1,1,10]:
w0 = np.zeros(2)
w,fiter,giter = newton_backtracking(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma),
lambda w: hF(w,X,I,J,gamma), w0, 1e-7, 1000)
plt.semilogy(giter)
iteration 1: f=41.588831, |g|=21.362260, step=1.000000 iteration 2: f=9.690587, |g|=5.491137, step=1.000000 iteration 3: f=3.890789, |g|=2.001470, step=1.000000 iteration 4: f=1.779552, |g|=0.737436, step=1.000000 iteration 5: f=1.015670, |g|=0.257076, step=1.000000 iteration 6: f=0.794580, |g|=0.073995, step=1.000000 iteration 7: f=0.760364, |g|=0.012390, step=1.000000 iteration 8: f=0.759077, |g|=0.000524, step=1.000000 iteration 9: f=0.759074, |g|=0.000001, step=1.000000 iteration 1: f=41.588831, |g|=21.362260, step=1.000000 iteration 2: f=10.147224, |g|=5.352513, step=1.000000 iteration 3: f=4.959125, |g|=1.802366, step=1.000000 iteration 4: f=3.565504, |g|=0.515521, step=1.000000 iteration 5: f=3.352376, |g|=0.085191, step=1.000000 iteration 6: f=3.344638, |g|=0.003410, step=1.000000 iteration 7: f=3.344625, |g|=0.000006, step=1.000000 iteration 1: f=41.588831, |g|=21.362260, step=1.000000 iteration 2: f=14.047299, |g|=4.186157, step=1.000000 iteration 3: f=12.023820, |g|=0.696488, step=1.000000 iteration 4: f=11.948721, |g|=0.029443, step=1.000000 iteration 5: f=11.948578, |g|=0.000058, step=1.000000 iteration 1: f=41.588831, |g|=21.362260, step=1.000000 iteration 2: f=28.834896, |g|=0.692093, step=1.000000 iteration 3: f=28.820387, |g|=0.002071, step=1.000000
On applique la méthode décrite dans la partie I (où l'on traitait des données en dimension $d=2$) à un jeu de données de dimension beaucoup plus grande. L'objectif est de distinguer des caractères manuscrits: pour cela, on a demandé à 6000 personnes d'écrire les chiffres de $0$ à $9$. Ici, l'on cherchera à distinguer les $0$ des $1$. Chaque point $x_i \in \Rsp^{784}$ est construit de la manière suivante:
La fonction read_mnist() ci-dessous télécharge ces données et retourne deux tableaux numpy images et labels:
# MNIST
import os
import struct
import urllib
import gzip
import shutil
try:
from urllib.request import urlretrieve
except ImportError:
from urllib import urlretrieve
def gunzip(iname,oname):
with gzip.open(iname, 'rb') as f_in:
with open(oname, 'wb') as f_out:
shutil.copyfileobj(f_in, f_out)
def download_mnist(dest):
src = 'http://yann.lecun.com/exdb/mnist'
for f in ['train-images-idx3-ubyte', 'train-labels-idx1-ubyte']:
urlretrieve('%s/%s.gz' % (src,f), '%s/%s.gz' % (dest,f))
gunzip('%s/%s.gz' % (dest,f),'%s/%s' % (dest,f))
def read_mnist():
d = os.path.expanduser('~/.m315-mnist-data')
if not os.path.isdir(d):
os.makedirs(d)
download_mnist(d)
with open('%s/train-labels-idx1-ubyte' % d, 'rb') as flbl:
magic, num = struct.unpack(">II", flbl.read(8))
labels = np.fromfile(flbl, dtype=np.int8)
with open('%s/train-images-idx3-ubyte' % d, 'rb') as fimg:
magic, num, rows, cols = struct.unpack(">IIII", fimg.read(16))
images = np.fromfile(fimg, dtype=np.uint8).reshape(len(labels), rows*cols)/255.
return images,labels
def show_mnist(image, cmap='gray'):
plt.figure()
imgplot = plt.imshow(np.reshape(image,(28,28)), interpolation='nearest', cmap=cmap)
images,labels = read_mnist()
Il faut passer quelques minutes à jouer avec les fonctions show_mnist() et print() (afficher les labels, et aussi les vecteurs images[i]) pour comprendre le jeu de données.
# exemple d'image: faire varier un peu i pour comprendre le jeu de données
i = 15
show_mnist(images[i])
print(labels[i])
7
QII.1 Construire un tableau $X0$ contenant $n=100$ images du chiffre $0$ extraite du tableau 'images', et un tableau $X1$ contenant également $n$ représentants du chiffre $1$ également extraits de images (X0 et X1 possèderont donc $n$ lignes et $784$ colonnes). Construire ensuite les données du problème d'optimisation:
ntot = len(labels)
n = 100 # nombre d'images par classe
X0 = [images[i] for i in range(ntot) if labels[i] == 0][0:n]
X1 = [images[i] for i in range(ntot) if labels[i] == 1][0:n]
X = np.vstack((X0,X1))
I = range(0,n)
J = range(n,2*n)
QII.2 Résoudre le problème d'optimisation $\min_{w\in\Rsp^d} F(w)$ pour $\gamma=0.1$ en utilisant la méthode de Newton.
gamma = .1
w0 = np.zeros(784)
w,fIter,gIter = newton_backtracking(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma),
lambda w: hF(w,X,I,J,gamma),w0,
err=1e-8,maxiter=5000)
iteration 1: f=138.629436, |g|=377.194291, step=1.000000 iteration 2: f=25.816351, |g|=89.958706, step=1.000000 iteration 3: f=9.159966, |g|=31.569106, step=1.000000 iteration 4: f=3.712922, |g|=11.593728, step=1.000000 iteration 5: f=1.701434, |g|=4.355502, step=1.000000 iteration 6: f=0.926405, |g|=1.640519, step=1.000000 iteration 7: f=0.640555, |g|=0.598944, step=1.000000 iteration 8: f=0.557049, |g|=0.190618, step=1.000000 iteration 9: f=0.543141, |g|=0.041410, step=1.000000 iteration 10: f=0.542359, |g|=0.003433, step=1.000000 iteration 11: f=0.542354, |g|=0.000029, step=1.000000
Soit $w$ la solution du problème d'optimisation et $u:x\mapsto\sigma(\sca{x}{w})$. L'idée de la régression logistique (cf intro du TD) est que $u(x) \simeq 1$ si $x$ représente un $1$ et $u(x) \simeq 0$ si $x$ représente un $0$. On n'a utilisé qu'une petite fraction des images pour construire $w$ (on a choisi $n=100$ images de $0$ alors que la base de données en contient $6000$): l'idée est de valider la méthode en comparant la classe ($0$ où $1$) devinée par la fonction $u$ à celle enregistrée dans le tableau labels sur l'ensemble des $0$ et $1$ de la base de données.
QII.3 Vérifier la pertinence de la méthode de la manière suivante :
total = 0
correct = 0
for i in range(ntot):
if labels[i] not in [0,1]:
continue
P = sigma(np.dot(images[i],w.T)) > .5
total = total+1
if P==labels[i]:
correct = correct+1
print('pourcentage de cas ou le chiffre 0 est correctement reconnu: %g%%' % (100.*float(correct)/float(total)))
pourcentage de cas ou le chiffre 0 est correctement reconnu: 99.6447%
QII.4: Afficher le vecteur $w$ sous la forme d'une image en utilisant show_mnist(). Interpréter.
show_mnist(np.reshape(w,(28,28)),cmap='RdBu')
QII.5 Adapter l'approche à la séparation des chiffres 0 de tout les autres chiffres. (on pourra prendre n=100 chiffres $0$ et $n$ chiffres de chaque autre classe, soit 1000 points de donnée au total).
n = 100 # nombre d'images par classe
X0 = [images[i] for i in range(ntot) if labels[i] == 0][0:n]
X1 = [images[i] for i in range(ntot) if labels[i] >= 1][0:9*n]
X = np.vstack((X0,X1))
I = range(0,n)
J = range(n,10*n)
gamma = 10
w0 = np.zeros(784)
w,fIter,gIter = newton_backtracking(lambda w: F(w,X,I,J,gamma),
lambda w: gF(w,X,I,J,gamma),
lambda w: hF(w,X,I,J,gamma),w0,
err=1e-8,maxiter=5000)
show_mnist(np.reshape(w,(28,28)),cmap='RdBu')
total = 0
correct = 0
for i in range(ntot):
R = labels[i] >= 1
P = sigma(np.dot(images[i],w.T)) > .5
total = total+1
if P==R:
correct = correct+1
print('pourcentage de cas ou le chiffre 0 est correctement reconnu: %g%%' % (100.*float(correct)/float(total)))
iteration 1: f=693.147181, |g|=2277.565327, step=1.000000 iteration 2: f=186.520113, |g|=572.343704, step=1.000000 iteration 3: f=98.207586, |g|=205.937616, step=1.000000 iteration 4: f=69.212164, |g|=71.872824, step=1.000000 iteration 5: f=61.692959, |g|=20.282336, step=1.000000 iteration 6: f=60.731338, |g|=3.026784, step=1.000000 iteration 7: f=60.706153, |g|=0.092316, step=1.000000 iteration 8: f=60.706130, |g|=0.000092, step=1.000000 pourcentage de cas ou le chiffre 0 est correctement reconnu: 98.6%