#!/usr/bin/env python # -*- coding: utf-8 -*- import numpy import scipy.spatial # データ読み込み f = open('nihonseiji.txt') head = f.readline() parties = head.strip().split('\t')[1:] vlist = [] for line in f: cols = line.strip().split('\t') matches = [float(col) for col in cols[1:]] vlist.append(matches) A = numpy.vstack(vlist) # 距離行列の計算 r, c = A.shape D = numpy.zeros((c, c)) VI = numpy.linalg.inv(numpy.cov(A)) for i in range(c): for j in range(i, c): u = A[:,i] v = A[:,j] # ユークリッド距離 D[i,j] = D[j,i] = scipy.spatial.distance.euclidean(u, v) # マハラノビス距離 # D[i,j] = D[j,i] = scipy.spatial.distance.mahalanobis(u, v, VI) # 相関行列 # print numpy.corrcoef(A.T)[i,j], # 分散共分散行列 # print numpy.cov(A) # データの個数 N = len(D) # 距離の2乗の行列 (arrayだと要素同士の掛け算になる) S = D * D # 中心化行列 one = numpy.eye(N) - numpy.ones((N,N))/N # ヤング・ハウスホルダー変換 P = - 1.0/2 * one * S * one # スペクトル分解 w,v = numpy.linalg.eig(P) ind = numpy.argsort(w) x1 = ind[-1] # 1番 x2 = ind[-2] # 2番 # print w[x1],w[x2] # 標準されたデータの固有値が求められているので標準偏差を掛けて座標を求める s = P.std(axis=0) w1 = s[x1] w2 = s[x2] X = [] Y = [] for i in range(N): X += [w1*v[i,x1]] Y += [w2*v[i,x2]] print parties[i], w1*v[i,x1], w2*v[i,x2]