nimfaでNMF
生命科学の分野でも、遺伝子発現情報やDNAメチル化情報のクラスター解析に使われているNMF(non-negative matrix factorization)を応用した手法を実装したPython用ライブラリにに、nimfaがあります。データマイニングソフトのOrangeにも採用されているようなので、使って見ました。
まず、こちらのGithubサイトから、レポジトリをzipファイルで取得して解凍。sudo python setup.py installで入れます。Python Image Library(PIL)が必要なので、なければこちらはpipなどで手軽に入れておきます。
NMFはやはりいくつのクラスに分類するかがもっとも効いてくるパラメータですので、これを2から9でふるサンプルスクリプトを載せておきます。といっても、estimate_rankというメソッドがあるので便利です。データは、行方向に変数、列方向にサンプルが並んだテキストファイルを、ほとんど神といっても過言では無いPandasのDataFrameとして読み込んでいます。一部、nimfaのサンプルコードを利用しています。

#!/usr/bin/env python

import pandas
import nimfa
import numpy as np
from scipy.cluster.hierarchy import linkage, leaves_list
from matplotlib.pyplot import savefig, imshow, set_cmap

def reorder(C):
"""
Reorder consensus matrix.

:param C: Consensus matrix.
:type C: `numpy.matrix`
http://nimfa.biolab.si/nimfa.examples.all_aml.html
"""
c_vec = np.array([C[i,j] for i in xrange(C.shape[0] - 1) for j in xrange(i + 1, C.shape[1])])
# convert similarities to distances
Y = 1 - c_vec
Z = linkage(Y, method = 'average')
# get node ids as they appear in the tree from left to right(corresponding to observation vector idx)
ivl = leaves_list(Z)
ivl = ivl[::-1]
return C[:, ivl][ivl, :],ivl

def plot(_C, file_name):
"""
Plot reordered consensus matrix.

:param C: Reordered consensus matrix.
:type C: `numpy.matrix`
:param rank: Factorization rank.
:type rank: `int`
http://nimfa.biolab.si/nimfa.examples.all_aml.html
"""
C,idx = reorder(_C)
imshow(np.array(C))
set_cmap("RdBu_r")
savefig(file_name)

def run_nmf(data):
V = np.matrix(data.as_matrix())
model = nimfa.mf(V, method = "lsnmf", max_iter = 200)
est_rank = model.estimate_rank(range=xrange(2,10),n_run=100)
for i in xrange(2,10):
print('rank = {}'.format(i))
print('cophenetic coef = {}'.format(est_rank[i]['cophenetic']))
print('sample prediction')
print(est_rank[i]['predict_samples'][0])
f_name = 'clustered_consensus_matrix_{}.png'.format(i)
plot(est_rank[i]['consensus'],f_name)
return est_rank

if __name__ == '__main__':
data = pandas.read_table('data.txt',index_col=0,na_values=(' '))
run_nmf(data)

Consensus行列を可視化するとこんな感じになります。
4クラスのConsensus行列を可視化
スポンサーサイト

テーマ:数学 - ジャンル:学問・文化・芸術

【2013/02/26 14:49】 | データ解析 | トラックバック(0) | コメント(0) | page top↑
<<RでROC | ホーム | numpyコンパイルエラー>>
コメント
コメントの投稿














管理者にだけ表示を許可する

トラックバック
トラックバックURL
→http://tanopy.blog79.fc2.com/tb.php/93-0776f04a
この記事にトラックバックする(FC2ブログユーザー)
| ホーム |