Valeurs manquantes et factorisation de matrices#

Links: notebook, html, PDF, python, slides, GitHub

Réflexion autour des valeur manquantes et de la factorisation de matrice positive.

from jyquickhelper import add_notebook_menu
add_notebook_menu()

%matplotlib inline

Matrice à coefficients aléatoires#

On étudie la factorisation d’une matrice à coefficients tout à fait aléatoires qui suivent une loi uniforme sur l’intervalle [0,1]. Essayons sur une petite matrice :

from numpy.random import rand
M = rand(3, 3)
M

array([[ 0.05119593,  0.43722929,  0.9290821 ],
       [ 0.4588466 ,  0.14187813,  0.23762633],
       [ 0.9768084 ,  0.47674026,  0.79044526]])

from sklearn.decomposition import NMF
mf = NMF(1)
mf.fit_transform(M)

array([[ 0.67825803],
       [ 0.38030919],
       [ 1.02295362]])

La matrice précédente est la matrice W dans le produit WH, la matrice qui suit est H.

mf.components_

array([[ 0.73190904,  0.50765757,  0.92611883]])

mf.reconstruction_err_ / (M.shape[0] * M.shape[1])

0.07236890712696428

On recalcule l’erreur :

d = M - mf.fit_transform(M) @ mf.components_
a = d.ravel()
e = a @ a.T
e ** 0.5 / (M.shape[0] * M.shape[1])

0.072368907126964283

e.ravel()

array([ 0.42421796])

Et maintenant sur une grande et plus nécessairement carrée :

M = rand(300, 10)
mf = NMF(1)
mf.fit_transform(M)
mf.reconstruction_err_ / (M.shape[0] * M.shape[1])

0.004996164872801101

L’erreur est la même :

errs = []
rangs = list(range(1, 11))
for k in rangs:
    mf = NMF(k)
    mf.fit_transform(M)
    e = mf.reconstruction_err_ / (M.shape[0] * M.shape[1])
    errs.append(e)

import pandas
df = pandas.DataFrame(dict(rang=rangs, erreur=errs))
df.plot(x="rang", y="erreur")

<matplotlib.axes._subplots.AxesSubplot at 0x199924d8668>
../_images/valeurs_manquantes_mf_16_1.png

Matrice avec des vecteurs colonnes corrélés#

Supposons maintenant que la matrice précédente M est de rang 3. Pour s’en assurer, on tire une matrice aléalatoire avec 3 vecteurs colonnes et on réplique des colonnes jusqu’à la dimension souhaitée.

from numpy import hstack
M = rand(300, 3)
M = hstack([M, M, M, M[:,:1]])
M.shape

(300, 10)

errs = []
rangs = list(range(1, 11))
for k in rangs:
    mf = NMF(k)
    mf.fit_transform(M)
    e = mf.reconstruction_err_ / (M.shape[0] * M.shape[1])
    errs.append(e)

import pandas
df = pandas.DataFrame(dict(rang=rangs, erreur=errs))
df.plot(x="rang", y="erreur")

<matplotlib.axes._subplots.AxesSubplot at 0x199925d6630>
../_images/valeurs_manquantes_mf_20_1.png

On essaye à nouveausur une matrice un peu plus petite.

M = rand(3, 2)
M = hstack([M, M[:,:1]])
M

array([[ 0.27190312,  0.6497563 ,  0.27190312],
       [ 0.44853292,  0.87097224,  0.44853292],
       [ 0.29424835,  0.65106952,  0.29424835]])

mf = NMF(2)
mf.fit_transform(M)

array([[ 0.61835197,  0.        ],
       [ 0.82887888,  0.29866219],
       [ 0.61960446,  0.07743224]])

mf.components_

array([[ 0.43972536,  1.05078419,  0.43972536],
       [ 0.28143493,  0.        ,  0.28143493]])

La dernière colonne est identique à la première.

Matrice identité#

Et maintenant si la matrice M est la matrice identité, que se passe-t-il ?

from numpy import identity
M = identity(3)
M

array([[ 1.,  0.,  0.],
       [ 0.,  1.,  0.],
       [ 0.,  0.,  1.]])

mf = NMF(1)
mf.fit_transform(M)

array([[ 0.],
       [ 1.],
       [ 0.]])

mf.components_

array([[ 0.,  1.,  0.]])

mf.reconstruction_err_ ** 2

2.0000000000000004

On essaye avec k=2.

mf = NMF(2)
mf.fit_transform(M)

array([[ 0.        ,  0.        ],
       [ 0.        ,  1.03940448],
       [ 0.95521772,  0.        ]])

mf.components_

array([[ 0.        ,  0.        ,  1.04688175],
       [ 0.        ,  0.96208937,  0.        ]])

mf.reconstruction_err_ ** 2

1.0

Avec des vecteurs normés et indépendants (formant donc une base de l’espace vectoriel), l’algorithme aboutit à une matrice W égale au k premiers vecteurs et oublie les autres.

Matrice identité et représentation spatiale#

Pour comprendre un peu mieux ce dernier exemple, il est utile de chercher d’autres solutions dont l’erreur est équivalente.

def erreur_mf(M, W, H):
    d = M - W @ H
    a = d.ravel()
    e = a @ a.T
    e ** 0.5 / (M.shape[0] * M.shape[1])
    return e

M = identity(3)
mf = NMF(2)
W = mf.fit_transform(M)
H = mf.components_
erreur_mf(M, W, H)

1.0

W

array([[ 0.        ,  0.        ],
       [ 0.9703523 ,  0.        ],
       [ 0.        ,  1.02721047]])

H

array([[ 0.        ,  1.03055354,  0.        ],
       [ 0.        ,  0.        ,  0.97351032]])

W @ H

array([[ 0.,  0.,  0.],
       [ 0.,  1.,  0.],
       [ 0.,  0.,  1.]])

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
wh = W @ H
ax.scatter(M[:,0], M[:,1], M[:,2], c='b', marker='o', s=20)
ax.scatter(wh[:,0], wh[:,1], wh[:,2], c='r', marker='^')

<mpl_toolkits.mplot3d.art3d.Path3DCollection at 0x19992d2a5c0>
../_images/valeurs_manquantes_mf_41_1.png

Et si on pose maintenant :

import numpy
W = numpy.array([[0.5, 0.5, 0], [0, 0, 1]]).T
H = numpy.array([[1, 1, 0], [0.0, 0.0, 1.0]])
W

array([[ 0.5,  0. ],
       [ 0.5,  0. ],
       [ 0. ,  1. ]])

H

array([[ 1.,  1.,  0.],
       [ 0.,  0.,  1.]])

W @ H

array([[ 0.5,  0.5,  0. ],
       [ 0.5,  0.5,  0. ],
       [ 0. ,  0. ,  1. ]])

erreur_mf(M, W, H)

1.0

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
wh = W @ H
ax.scatter(M[:,0], M[:,1], M[:,2], c='b', marker='o', s=20)
ax.scatter(wh[:,0], wh[:,1], wh[:,2], c='r', marker='^')

<mpl_toolkits.mplot3d.art3d.Path3DCollection at 0x19992a2e5f8>
../_images/valeurs_manquantes_mf_47_1.png

On peut voir la matrice M comme un ensemble de n points dans un espace vectoriel. La matrice W est un ensemble de k < n points dans le même espace. La matrice WH, de rang k est une approximation de cet ensemble dans le même espace, c’est aussi n combinaisons linéaires de k points de façon à former n points les plus proches proches de n points de la matrice M.