scipy を用いた Latent semantic indexing

#coding: utf-8

import numpy as np
from scipy import linalg

'''
Latent semantic indexing
'''
def main():
    # P.Baldi et al., 確率モデルによるWebデータ解析法, 森北出版, pp.96-98 に
    # 掲載されている行列を用いる
    X = np.matrix([[1,1,0,0,0,0,0,0,0,0,0],
                   [0,0,0,1,1,0,0,0,0,0,0],
                   [0,0,0,1,1,1,0,0,0,0,0],
                   [0,1,1,1,0,0,0,0,0,0,0],
                   [1,0,0,0,0,1,1,0,0,0,0],
                   [0,0,0,0,0,0,0,1,1,0,0],
                   [0,0,0,0,0,0,0,0,0,1,2],
                   [0,0,0,0,0,0,1,0,0,1,0],
                   [0,0,0,0,0,0,0,0,1,1,0],
                   [0,0,0,0,0,0,0,1,0,0,1]])

    # scipy.linalg.svd
    # http://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.svd.html
    u, sigma, v = linalg.svd(X)

    rank = np.shape(X)[0]
    u = np.matrix(u)
    # linalg.svd() の結果,sigma には特異値のリスト(降順)が返ってくるため
    # linalg.diagsvd を用いて行列に変換する
    # http://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.diagsvd.html
    sigma = np.matrix(linalg.diagsvd(sigma, rank, rank))
    v = np.matrix(v)

    # 特異値の大きい方から2個をとって X を再構成
    z = 2
    u2 = u[:, :z]
    sigma2 = sigma[:z, :z]
    v2 = v[:z, :]
    Xhat = u2 * sigma2 * v2

    print Xhat

    np.savetxt('result.txt', Xhat, fmt='%.02f')
    
if __name__ == '__main__':
    main()

実行結果は以下のとおり:

[[  1.18742376e-01   1.48468092e-01   1.02079902e-01   3.36802551e-01
    2.34722649e-01   2.07746172e-01   9.14487067e-02  -2.08551485e-03
    3.41743214e-03   2.40646404e-02   6.94193966e-04]
 [  2.52738635e-01   3.18786565e-01   2.19455902e-01   7.23753688e-01
    5.04297787e-01   4.44151468e-01   1.87007032e-01  -2.02395303e-02
   -4.72794523e-03   1.21686493e-02  -5.85747967e-02]
 [  3.45097561e-01   4.34173811e-01   2.98781826e-01   9.85492384e-01
    6.86710559e-01   6.05672638e-01   2.58391353e-01  -2.13348466e-02
   -1.66967515e-03   3.21884065e-02  -5.60326060e-02]
 [  2.59749943e-01   3.27600603e-01   2.25520701e-01   7.43758428e-01
    5.18237727e-01   4.56451866e-01   1.92276046e-01  -2.06330486e-02
   -4.73153071e-03   1.29213297e-02  -5.95614097e-02]
 [  1.87913028e-01   2.30024227e-01   1.57670040e-01   5.20783588e-01
    3.63113547e-01   3.25264417e-01   1.58277142e-01   2.47446747e-02
    2.67106562e-02   1.07397837e-01   1.07687758e-01]
 [  6.67965446e-03  -5.34773717e-03  -5.02206960e-03  -1.49947726e-02
   -9.97270298e-03   1.96295375e-03   4.28127228e-02   7.78060127e-02
    5.93812774e-02   1.93945477e-01   2.96197591e-01]
 [  4.62012230e-02  -2.07481947e-02  -2.19752154e-02  -6.34780797e-02
   -4.15028644e-02   2.51041386e-02   2.51467991e-01   4.45785034e-01
    3.40555624e-01   1.11315083e+00   1.69762265e+00]
 [  6.18128209e-02   5.37005262e-02   3.46060061e-02   1.16890844e-01
    8.22848376e-02   9.14040781e-02   1.12458081e-01   1.33074100e-01
    1.03684099e-01   3.44110267e-01   5.10254277e-01]
 [  2.34484339e-02   4.03363860e-03   2.90538891e-04   3.86562148e-03
    3.57508259e-03   2.30780273e-02   8.75820325e-02   1.43408667e-01
    1.09918087e-01   3.60212333e-01   5.46747060e-01]
 [  1.59795054e-02  -1.73708263e-02  -1.56109790e-02  -4.72126530e-02
   -3.16016740e-02   1.44687434e-03   1.15006052e-01   2.12170201e-01
    1.61833403e-01   5.28322325e-01   8.07542678e-01]]

読みづらいから,小数点以下2桁で表示すると次のようになる。

0.12  0.15  0.10  0.34  0.23  0.21  0.09 -0.00  0.00  0.02  0.00
0.25  0.32  0.22  0.72  0.50  0.44  0.19 -0.02 -0.00  0.01 -0.06
0.35  0.43  0.30  0.99  0.69  0.61  0.26 -0.02 -0.00  0.03 -0.06
0.26  0.33  0.23  0.74  0.52  0.46  0.19 -0.02 -0.00  0.01 -0.06
0.19  0.23  0.16  0.52  0.36  0.33  0.16  0.02  0.03  0.11  0.11
0.01 -0.01 -0.01 -0.01 -0.01  0.00  0.04  0.08  0.06  0.19  0.30
0.05 -0.02 -0.02 -0.06 -0.04  0.03  0.25  0.45  0.34  1.11  1.70
0.06  0.05  0.03  0.12  0.08  0.09  0.11  0.13  0.10  0.34  0.51
0.02  0.00  0.00  0.00  0.00  0.02  0.09  0.14  0.11  0.36  0.55
0.02 -0.02 -0.02 -0.05 -0.03  0.00  0.12  0.21  0.16  0.53  0.81