[[ノート/ノート]]~
[[ノート/R]]~
訪問者数 &counter();      最終更新 &lastmod();

**Rで確率モデルによるクラスタリング(2010-07-21)  [#w48d43a8]
教科書: [[新納浩幸:Rで学ぶクラスタ解析(オーム社、ISBN978-4-274-06703-7):http://www.amazon.co.jp/R%E3%81%A7%E5%AD%A6%E3%81%B6%E3%82%AF%E3%83%A9%E3%82%B9%E3%82%BF%E8%A7%A3%E6%9E%90-%E6%96%B0%E7%B4%8D-%E6%B5%A9%E5%B9%B8/dp/4274067033/ref=sr_1_1?ie=UTF8&s=books&qid=1277623068&sr=8-1]]の第6章~
   [[>金明哲先生:「Rによるデータサイエンス - データ解析の基礎から最新手法まで」:http://www.morikita.co.jp/shoshi/ISBN978-4-627-09601-1.html]]  [[>アマゾンのページ:http://www.amazon.co.jp/R%E3%81%AB%E3%82%88%E3%82%8B%E3%83%87%E3%83%BC%E3%82%BF%E3%82%B5%E3%82%A4%E3%82%A8%E3%83%B3%E3%82%B9-%E3%83%87%E3%83%BC%E3%82%BF%E8%A7%A3%E6%9E%90%E3%81%AE%E5%9F%BA%E7%A4%8E%E3%81%8B%E3%82%89%E6%9C%80%E6%96%B0%E6%89%8B%E6%B3%95%E3%81%BE%E3%81%A7-%E9%87%91-%E6%98%8E%E5%93%B2/dp/4627096011]]~
   [[>金明哲先生の教材ページ:http://mjin.doshisha.ac.jp/R/29/29.html]]  [[>シリーズのトップページ:http://mjin.doshisha.ac.jp/R/]]~
   [[>データマイニング+WEB勉強会のページ:http://d.hatena.ne.jp/hamadakoichi/20100415/p2]]

   [[>R用パッケージMclustのページ(Washington大 統計学部) (http://www.stat.washington.edu/mclust/):http://www.stat.washington.edu/mclust/]]

   [[>RではないGMM (Gaussian Mixture Model) のパッケージ:http://mlpack.org/doxygen.php?doc=namespacemlpack_1_1gmm.html]] 〜 mlpack(machine learning package, C++のパッケージ)中のgmmライブラリ
   [[>RではないGMM (Gaussian Mixture Model) のパッケージ:http://mlpack.org/doxygen.php?doc=namespacemlpack_1_1gmm.html]] 〜 mlpack(machine learning package, C++のパッケージ)中のgmmライブラリ~
     [[Apache Mahout:http://mahout.apache.org//]]とその[[アルゴリズム:https://cwiki.apache.org/confluence/display/MAHOUT/Algorithms]]~
     [[Shogun toolbox:http://www.shogun-toolbox.org/]]~


***理論編(混合分布とそのEMアルゴリズムによる推定) [#p71cb231]
-[[混合分布:http://ibisforest.org/index.php?%E6%B7%B7%E5%90%88%E5%88%86%E5%B8%83]]
--教科書(1)「データマイニングの基礎」(元田浩etal) オーム社 3.2の6 (p85)
--教科書(2)「パターン認識と機械学習」(C.M.ビショップ、元田浩etal訳) Springer 第9章

-[[金田尚久・新居玄武:混合分布問題 その基礎からカーネル降下法まで:http://www.gakushuin.ac.jp/univ/eco/gakkai/pdf_files/keizai_ronsyuu/contents/4601/4601arai/4601arai.htm]]

-[[正田備也:ガウス混合分布のEMアルゴリズムによるパラメータ推定:http://www.iris.dti.ne.jp/~tmasada/2007032501.pdf]]

-[[田中洸一:混合分布モデル:http://nlp.dse.ibaraki.ac.jp/~shinnou/zemi2008/Rclustering/r-tanaka-0617.pdf]]

-[[江口真透:混合分布:http://www.ism.ac.jp/~eguchi/pdf/openmix2000.pdf]]

-[[田中研太郎:混合分布モデルにおける一致推定量の構成:http://www.i.u-tokyo.ac.jp/coe/report/H15/21COE-ISTSC-H15_4_5_11.pdf]]

-[[EMアルゴリズムを用いた混合分布の推定:http://mist.murase.m.is.nagoya-u.ac.jp/document/group__mixture__group.html]]

-[[赤穂昭太郎 ECM 法を用いた確率分布の位置, 尺度, 回転パラメータの推定法:http://staff.aist.go.jp/s.akaho/papers/IEICE99.pdf]]

-(計算パッケージ)
--[[混合分布に対する EM アルゴリズム (テスト バージョン):http://staff.aist.go.jp/s.akaho/index-j.html]]
--[[OpenCVにあるEMアルゴリズム:EMアルゴリズム(Expectation-Maximization):http://opencv.jp/opencv-2svn/cpp/expectation-maximization.html]]

--[[sage/PRML - 混合ガウス分布のEMアルゴリズム:http://www.pwv.co.jp/~take/TakeWiki/index.php?sage%2FPRML%20-%20%E6%B7%B7%E5%90%88%E3%82%AC%E3%82%A6%E3%82%B9%E5%88%86%E5%B8%83%E3%81%AEEM%E3%82%A2%E3%83%AB%E3%82%B4%E3%83%AA%E3%82%BA%E3%83%A0]]

-[[Tutorial of the EM Algorithm and its Application to Parameter Estimation for Gaussian Mixture and Hidden Markov Models:http://lasa.epfl.ch/teaching/lectures/ML_Phd/Notes/GP-GMM.pdf]]

-[[Dr. Dowe's page on Mixture Modelling:http://www.csse.monash.edu.au/~dld/mixture.modelling.page.html]]

-[[ EMアルゴリズムの改良:http://hillbig.cocolog-nifty.com/do/2009/04/em-cc07.html]]、 [[http://d.hatena.ne.jp/isseing333/20110412/1302619320]]

-[[名古屋大能登原:http://www.nsc.nagoya-cu.ac.jp/~noto/emalgo.pdf]]

-[[北大井上:http://chaosweb.complex.eng.hokudai.ac.jp/~j_inoue/NC_BN2004OCT_inoue.pdf]]

-[[高知大本田:http://www.is.kochi-u.ac.jp/~honda/studentarc/03sahara.pdf]]

-[[茨城大新納:http://nlp.dse.ibaraki.ac.jp/~shinnou/zemi2007/EM/EM-abe.pdf]]

-[[EMアルゴリズムの教科書:http://d.hatena.ne.jp/isseing333/20110412/1302619320]]

-[[機械学習:http://www40.atwiki.jp/geiinbashoku/pages/22.html]]

-[[多変量混合正規分布のEMアルゴリズム・サンプルプログラム その1:http://ameblo.jp/p630/entry-10962467573.html]]


**Rによる解析5 (山の具合を描いてみた) (2010/07/21 17:30) [#da1fd9b1]
山の形に見えるのが結構難しい。~
3次元プロットすると、X-Yの解像度が高いと、細い縦の線になってしまう。つまり表面が平滑化されていない。そこで、X-Yの解像度を256, 128, 64で比べてみる。

[A] 少ないデータの方(FSC-SSC_100113_24h_40C_+20-1)~
|CENTER: X-Y軸解像度=64|CENTER: X-Y軸解像度=128|
|&ref(HillShape_64_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|&ref(HillShape_128_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|

[B] 多いデータの方(FSC-SSC_100113_24h_40C_+9-1)~
|CENTER: X-Y軸解像度=64|CENTER: X-Y軸解像度=128|
|&ref(HillShape_64_FSC-SSC_100113_24h_40C_+9-1.png,,400x400);|&ref(HillShape_128_FSC-SSC_100113_24h_40C_+9-1.png,,400x400);|

[C]自動クラスタ化の結果との比較~
上記を踏まえて、X-Y解像度を64に設定した図で考えることにする。

少ないデータの方(FSC-SSC_100113_24h_40C_+20-1)しか自動クラスタ化をしていないので、それで比較する。~
|CENTER: X-Y軸解像度=64|CENTER: 自動クラスタ化の結果(クラスタ数=6)|
|&ref(HillShape_64_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|&ref(classification_20000_6_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|
左側の図から山を読み取ると、
-X=36付近に高さ80程度のピーク、Xを1-1024のスケールに直すと16倍してX=576。これは右図の緑に対応するだろう。
-X=38付近に高さ80程度のピーク、Xを1-1024のスケールに直すと16倍してX=608。これは右図の青に対応するだろう。
-X=42〜49付近に高さ30程度のピーク、Xを1-1024のスケールに直すと16倍してX=672〜784。これは右図の赤に対応するだろう。


**Rによる解析4 (FL1もクラスタ化してみた) (2010/07/16 18:30) [#g837fc3a]
FSC-FL1は、自動クラスタ化で、クラスタ数は6になった。

|                              |CENTER: FSC-FL1 |CENTER: FSC-SSC (6クラスタ)  |
|総データ数                    |CENTER: 20,000  |CENTER: 20,000  |
|ノイズフィルタ処理後のデータ数|CENTER:  3,245  |CENTER: 4,562   |
| クラスタの数                 |CENTER:  6      |CENTER:  6      |
||&ref(classification_20000_FSC-FL1_100113_24h_40C_+20-1.png,,400x400);|&ref(classification_20000_6_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|
   
| 総データ数        |>|>|>|>|>|CENTER:FSC-FL1|              |>|>|>|>|>|CENTER:FSC-SSC|
|  クラスタ番号      |1  |2  |3  |4  |5   |6  |    |1  |2  |3  |4   |5  |6  |
|(色別)        |青|赤|黄緑|空色|ピンク|濃緑| |青|赤|黄緑|空色|ピンク|濃緑|
|                |>|>|>|>|>|>|>|>|>|>|>|>|CENTER:両図のクラスタ間の対応は未決定なので、要注意|
| クラスタの大きさ(割合)|0.23|0.27|0.30|0.42|0.03|0.13| |0.37|0.29|0.21|0.02|0.02|0.08|
| 同上  (点数)| 752| 876| 973|135| 92|417| |1694|1308| 951| 102|110|398|
| クラスタの平均 x| 712| 617| 567|998|408|904| | 602| 690| 545| 999|376|904|
| 同上 y          | 631| 531| 514|978|289|981| | 586| 666| 560|1012|253|941|
| クラスタの分散 x-x|2129| 951|2375|164| 9637|18| |947|3150|2791|190| 6395|14|
| 同上 y-x          |1654| 983|1046| -1| -274|-1| |321|2193| 937|106|-1959| 4|
| 同上 y-y          |1841|2411|1598|  1|15453| 2| |385|1983| 990| 63| 5303| 7|

FSC-FL1の出力データの詳細
 $pro
 [1] 0.23167330 0.26984454 0.29997641 0.04160247 0.02839789 0.12850539
 $mean
          [,1]     [,2]     [,3]     [,4]     [,5]     [,6]
 [1,] 712.3932 616.5135 567.3878 998.4148 408.2780 903.6547
 [2,] 631.3966 531.1113 514.4008 978.0889 288.7181 981.1727
 $variance
 $variance$modelName
 [1] "VVV"
 $variance$d
 [1] 2
 $variance$G
 [1] 6
 $variance$sigma
 1
          [,1]     [,2]
 [1,] 2129.354 1654.443
 [2,] 1654.443 1840.998
 2
          [,1]      [,2]
 [1,] 950.7647  982.5923
 [2,] 982.5923 2411.1313
 3
          [,1]     [,2]
 [1,] 2374.617 1046.110
 [2,] 1046.110 1597.614
 4
            [,1]       [,2]
 [1,] 163.813114 -1.2368724
 [2,]  -1.236872  0.4661728
 5
           [,1]       [,2]
 [1,] 9637.2769  -273.5086
 [2,] -273.5086 15453.0124
 6
            [,1]       [,2]
 [1,] 17.5594086 -0.6310232
 [2,] -0.6310232  2.3634733
 
 [1] "Total data counts:  6490"
 [1] "Cluster sizes:"
 [1] 1504 1751 1947  270  184  834
          [,1]     [,2]     [,3]     [,4]     [,5]     [,6]
 [1,] 712.3932 616.5135 567.3878 998.4148 408.2780 903.6547
 [2,] 631.3966 531.1113 514.4008 978.0889 288.7181 981.1727



**Rによる解析3(点数を増やしてみると) (2010/07/16) [#ob7c3581]

***[A] クラスタ数自動決定(=8)の場合 [#t70aeac3]
計算時間を余分にかけて、測定点数を増やしてみる。同一実験データ(100113 24h 40C +20-1_nonlimit.csv)で、10,000点から20,000点へ。クラスタ数自動決定をすると8になった。
|総データ数                    |CENTER: 20,000  |CENTER: 10,000  |
|ノイズフィルタ処理後のデータ数|CENTER: 4,562   |CENTER: 2,592   |
| クラスタの数                 |CENTER:  8      |CENTER:  6      |
||&ref(classification_20000_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|&ref(classification_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|
   

| 総データ数        |>|>|>|>|>|>|>|CENTER:20,000|              |>|>|>|>|>|CENTER:10,000|
|  クラスタ番号      |1  |2  |3  |4  |5   |6  |7  |8  |    |1  |2  |3  |4   |5  |6  |
|(色別)        |青|赤|黄緑|空色|ピンク|濃緑|紫|黄土| |青|赤|黄緑|空色|ピンク|濃緑|
| クラスタの大きさ(割合)|0.11|0.31|0.27|0.11|0.02|0.02|0.09|0.07||0.36|0.21|0.24|0.04|0.03|0.12|
| 同上  (点数)|1028|2823|2421|1014| 204| 217| 796| 621|| 942| 538| 622| 92| 76| 322|
| クラスタの平均 x|619|586|659|524|999 |377|903|748||608|547|699|1000|357|904|
| 同上 y        |580|585|643|546|1012|249|941|700||589|568|674|1013|255|941|
| クラスタの分散 x-x|526|1234|2978|3026|190|6398| 14|966||582|3426|2525|189|7504|13|
| 同上 y-x          |468| 361|2377| 775|106|-1724| 4|899||250|1041|1743|106|-1003| 4|
| 同上 y-y          |537| 345|2438| 875| 63|4444|  7|989||353| 931|1677| 64|2800| 6|

***[B] クラスタ数固定(=6)の場合 [#tebe6e6c]
一方で、クラスタ数を10,000点時と同じ6つに固定すると、次のようなクラスタを作る。

|総データ数                    |CENTER: 20,000  |CENTER: 10,000  |
|ノイズフィルタ処理後のデータ数|CENTER: 4,562   |CENTER: 2,592   |
| クラスタの数                 |CENTER:  6      |CENTER:  6      |
||&ref(classification_20000_6_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|&ref(classification_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);|
   

| 総データ数        |>|>|>|>|>|CENTER:20,000|              |>|>|>|>|>|CENTER:10,000|
|  クラスタ番号      |1  |2  |3  |4  |5   |6  |    |1  |2  |3  |4   |5  |6  |
|(色別)        |青|赤|黄緑|空色|ピンク|濃緑| |青|赤|黄緑|空色|ピンク|濃緑|
|                |>|>|>|>|>|>|>|>|>|>|>|>|CENTER:クラスタ2(赤)と3(黄緑)が逆転しているので注意|
| クラスタの大きさ(割合)|0.37|0.29|0.21|0.02|0.02|0.08||0.36|0.21|0.24|0.04|0.03|0.12|
| 同上  (点数)|1694|1308| 951| 102|110|398|| 942| 538| 622| 92| 76| 322|
| クラスタの平均 x|602|690|545| 999|376|904||608|547|699|1000|357|904|
| 同上 y          |586|666|560|1012|253|941||589|568|674|1013|255|941|
| クラスタの分散 x-x|947|3150|2791|190| 6395|14||582|3426|2525|189|7504|13|
| 同上 y-x          |321|2193| 937|106|-1959| 4||250|1041|1743|106|-1003| 4|
| 同上 y-y          |385|1983| 990| 63| 5303| 7||353| 931|1677| 64|2800| 6|


***[A] クラスタ数自動決定の場合の出力データ [#hd4a348b]
20,000点時、クラスタ数=8(自動決定)の出力データ
 $pro
 [1] 0.11269930 0.30937046 0.26537431 0.11109090 0.02235861 0.02375630
 [7] 0.08724244 0.06810768
 $mean
          [,1]     [,2]     [,3]     [,4]      [,5]     [,6]     [,7]
 [1,] 619.0596 585.5307 658.6493 523.5670  998.6863 376.6446 903.6633
 [2,] 579.8497 584.8812 643.2790 546.4888 1011.5196 248.8440 940.8317
          [,8]
 [1,] 748.1149
 [2,] 700.2097
 $variance$modelName
 [1] "VVV"
 $variance$d
 [1] 2
 $variance$G
 [1] 8
 $variance$sigma
 1
          [,1]     [,2]
 [1,] 526.1028 467.6122
 [2,] 467.6122 537.2284
 2
           [,1]     [,2]
 [1,] 1234.2949 361.0131
 [2,]  361.0131 345.3007
 3
          [,1]     [,2]
 [1,] 2977.924 2376.559
 [2,] 2376.559 2437.739
 4
          [,1]     [,2]
 [1,] 3026.429 774.5300
 [2,]  774.530 874.8733
 5
          [,1]     [,2]
 [1,] 190.1761 105.9473
 [2,] 105.9473  63.2104
 6
           [,1]      [,2]
 [1,]  6397.935 -1723.837
 [2,] -1723.837  4444.096
 7
           [,1]     [,2]
 [1,] 14.062524 3.935784
 [2,]  3.935784 6.748043
 8
          [,1]     [,2]
 [1,] 965.6994 899.0406
 [2,] 899.0406 988.8561
 
 [1] "Total data counts:  9124"
 [1] "Cluster sizes:"
 [1] 1028 2823 2421 1014  204  217  796  621
          [,1]     [,2]     [,3]     [,4]      [,5]     [,6]     [,7]
 [1,] 619.0596 585.5307 658.6493 523.5670  998.6863 376.6446 903.6633
 [2,] 579.8497 584.8812 643.2790 546.4888 1011.5196 248.8440 940.8317
          [,8]
 [1,] 748.1149
 [2,] 700.2097

***[B] クラスタ数を6に固定した場合の出力データ [#i9faf7d5]
クラスタ数を10,000点時と同じ6つに固定したときのパラメタ出力。
 $pro
 [1] 0.37127726 0.28670306 0.20830639 0.02235861 0.02411224 0.08724244
 $mean
          [,1]     [,2]     [,3]      [,4]     [,5]     [,6]
 [1,] 601.7230 689.6276 545.2734  998.6863 375.6309 903.6633
 [2,] 585.6392 666.4851 560.2669 1011.5196 252.5498 940.8317
 $variance
 $variance$modelName
 [1] "VVV"
 $variance$d
 [1] 2
 $variance$G
 [1] 6
 
 $variance$sigma
 1
          [,1]     [,2]
 [1,] 947.0633 321.2171
 [2,] 321.2171 385.4888
 2
          [,1]     [,2]
 [1,] 3149.918 2193.402
 [2,] 2193.402 1983.022
 3
          [,1]     [,2]
 [1,] 2790.827 936.6290
 [2,]  936.629 990.1658
 4
          [,1]     [,2]
 [1,] 190.1761 105.9473
 [2,] 105.9473  63.2104
 5
           [,1]      [,2]
 [1,]  6394.872 -1959.228
 [2,] -1959.228  5302.817
 6
           [,1]     [,2]
 [1,] 14.062524 3.935784
 [2,]  3.935784 6.748043
 
 [1] "Total data counts:  9124"
 [1] "Cluster sizes:"

 [1] 1694 1308  951  102  110  398
          [,1]     [,2]     [,3]      [,4]     [,5]     [,6]
 [1,] 601.7230 689.6276 545.2734  998.6863 375.6309 903.6633
 [2,] 585.6392 666.4851 560.2669 1011.5196 252.5498 940.8317


**Rによる解析2(雑音除去の効果と、Mclustでの表示法の改善) (2010/07/14) [#b9125016]
雑音除去がかなり効果を挙げることがわかった。具体的には、(1)プロットチャートの左端の線(Xの値が0)と、(2)全体に散らばっている≪おそらくランダムノイズと思われる点≫、の2種類を除去すると、余分なヤマが出なくなる。

例:100113 24h 40C +20-1_nonlimit.csvを入力とする場合
-全体で20,000個のデータがあるが、時間短縮のため最初の10,000個のデータのみ利用した
-それに対して、上記の(1),(2)での除去を行うと結果は2,592個のデータが残った。かなり少なくなったと言えるだろう。
-その2,592個のデータに対して、Mclustを用いて(正規分布として)ヤマを推定した。
-結果は
 最適のモデルはVVV、ヤマの数は6つ。
 
 それぞれの分布(ヤマ)の割合
 $pro
 [1] 0.36346510 0.20738490 0.24009412 0.03549383 0.02933365 0.12422840
 
 それぞれの分布(ヤマ)のデータ個数 (2,592×$pro)
 [1] 942 538 622  92  76 322
 
 それぞれの分布(ヤマ)の中心(平均値)
          [,1]     [,2]    [,3]      [,4]     [,5]    [,6]
 [1,] 608.3974 547.3208 699.312  999.8478 356.7617 903.913
 [2,] 589.3498 567.8984 673.731 1012.5870 254.5305 941.149
 
 それぞれの分布(ヤマ)の分散(sigma)
 1                        2                        3
         [,1]     [,2]             [,1]     [,2]            [,1]     [,2]
 [1,] 581.7725 250.0036   [1,] 3425.921 1040.665   [1,] 2524.561 1742.984
 [2,] 250.0036 353.2514   [2,] 1040.665  931.292   [2,] 1742.984 1677.382
 
 4                        5                         6
          [,1]      [,2]            [,1]      [,2]            [,1]     [,2]
 [1,] 188.6073 106.19802  [1,]  7504.097 -1002.931  [1,] 12.849581 4.329733
 [2,] 106.1980  63.98157  [2,] -1002.931  2800.267  [2,]  4.329733 6.114425

また、表示法として、
 mclust2Dplot(data=z1, what="classification", identify=TRUE, parameters=x$paramters, z=x$z, scale=TRUE)
のような関数を使うことが出来る。結果は、

&ref(classification_FSC-SSC_100113_24h_40C_+20-1.png,,400x400);

色の「凡例」が表示されていないのでわかりづらいが、~
  グループ1=青、2=赤、3=黄緑、4=空色、、5=ピンク、6=濃い緑~
に対応する。

真ん中の赤(グループ2、データ点数538)と青(グループ1、データ点数942)がメイン、
黄緑(グループ3、データ点数622)が右上に新しく伸びた別の山、
右上の濃い緑(グループ6、データ点数322)が定点用ビーズ、
更に右上の空色(グループ4、データ点数92)と左下のピンク(グループ5、データ点数76)
はたぶんノイズと思われる。

このように、かなりきれいな形で分布と測定点の帰属が推定できると思われる。他のデータについてもテストしてみたい。

**ゼミで使った資料  (2010/07/12追記) [#lccedfa8]
|&ref(スライド21.PNG,,400x300);|&ref(スライド22.PNG,,400x300);|


***Rによる解析(新納の6.4節)(2010/06/29) [#c3d85181]
手順
-パッケージmclustをインストール~
  パッケージ → パッケージのインストール → mclustを選択
-プログラム例は
 library(mclust)
 x <- matrix(scan("myiris.dat"), ncol=4, byrow=TRUE)
 mc <- Mclust(x, G=3,modelNames="EII")
 ans <- mc$classification
 pa <- mc$parameters
 pa$pro
 pa$mean
 plot(Mclust(x, G=1:20,modelNames="VEI"))  # modelNames="EEI","VII", etc
ここで、myiris.datはR中のサンプルデータiris= (花びらの長さ・幅・がく片の長さ・幅・irisの種類)が150データ=からはじめの4値を抜き出したもの。
 5.1 3.5 1.4 0.2
 4.9 3.0 1.4 0.2
 4.7 3.2 1.3 0.2
 4.6 3.1 1.5 0.2
 5.0 3.6 1.4 0.2
 ...

***フローサイトメトリー(大腸菌)データでのモデルベースクラスタリング [#w1c96f86]
実験1 〜 測定データ FSC-SSCをそのまま入力してモデルベースのクラスタリングをしてみる。

|&ref(スライド23.PNG,,400x300);|

データはほとんどそのまま使える。
 install.packages("mclust")
 library(mclust)
 # Read the data
 ##u <- matrix(scan("100113 24h 40C +9-1_nonlimit.csv", sep=","), ncol=8, byrow=TRUE)
 u <- matrix(scan("100113 24h 40C +20-1_nonlimit.csv", sep=","), ncol=8, byrow=TRUE)
   # ncol should be 8 because lines end with comma (with one more hidden column.)
 u <- u[,1:7]   # erase the hidden=empty column
 
mclustの使い方については、いろいろ書いてあるところもあるが、全体を見るにはオリジナルのドキュメントも良い。--> [[tr504.pdf:http://www.stat.washington.edu/research/reports/2006/tr504.pdf]]

 # The document starts with 
 BIC <- mclustBIC(u[1:6000,1:2])

uのうちの1:2 (FSCとSSC)を使ってみる。また、全データ(10000点)を使うとメモリ不足になるので、最初の4000点を使ってみる。Windowsで6000点ぐらいまでは動いた。なお、メモリスペースについては
 help(memory.size)
で説明が出てくるが、memory.size(NA) と打つと現在のメモリサイズ(Windows XPの場合デフォルトで2GB、最大で3GB、Win7では4GBらしい)が表示され、
 memory.limit(3071)
とすると上限を3GBに設定することが出来る。いずれにせよ、この程度では、全ポイントを含めると、メモリが足りなくなる。

結果で得られたBICを印刷すると以下の通りで、プロットすると図の通り。
  BIC:
         EII       VII       EEI       VEI       EVI       VVI       EEE
 1 -162779.0 -162779.0 -162262.9 -162262.9 -162262.9 -162262.9 -155018.1
 2 -156268.7 -155705.2 -155662.7 -155286.2 -155283.2 -155294.3 -153774.0
 3 -154328.1 -153104.1 -154337.0 -153033.5 -153878.3        NA -149468.2
 4 -151256.3 -149930.7 -151269.8 -149849.3 -150173.3        NA -148455.5
 5 -149347.0 -147278.0 -149355.8 -147088.4 -147818.7        NA -148300.4
 6 -148476.1 -145945.0 -148481.6 -145875.4 -146298.5        NA -147772.1
 7 -148497.6 -144866.6 -148494.2 -144805.9 -145657.5        NA -147800.0
 8 -148454.1 -144143.6 -148467.9 -144108.7 -145461.2        NA -147827.8
 9 -146690.5 -143815.6 -146592.5 -143822.5 -144603.3        NA -147842.7
         EEV       VEV       VVV
 1 -155018.1 -155018.1 -155018.1
 2 -150277.4 -148559.5 -148223.5
 3 -147289.1 -145205.7        NA
 4 -145856.6 -144873.9        NA
 5 -144800.7 -143775.9        NA
 6 -144804.5 -142855.7        NA
 7 -144002.2 -142175.9        NA
 8 -143779.0 -142215.5        NA
 9 -143724.4 -141796.0        NA
&ref(BIC_FSC-SSCplot.png);~
要するに、方法はVEVが最もよくて、クラスター数は9ということである。

|&ref(スライド24.PNG,,400x300);|&ref(スライド25.PNG,,400x300);|

出力した結果のオリジナル表記は下記の通りになっていた。
 > summary(BIC, data=u[1:6000,1:2],modelNames="VEV")
 classification table:
    1    2    3    4    5    6    7    8    9 
  117 1231  639  431 1013 1400  443  591  135 
 
 best BIC values:
     VEV,9     VEV,7     VEV,8 
 -141796.0 -142175.9 -142215.5 
となっていた。再度
 x <- Mclust(u[1:6000,1:2], G=9, modelNames="VEV")
を処理して、推定したパラメタを見てみると、まずそれぞれのクラスタの大きさαiは
 $pro
 [1] 0.06943546 0.20443535 0.09925967 0.09049414 0.13425212 0.17498290
 [7] 0.10774480 0.09694017 0.02245540
で、i=2, 3, 5, 6, 8などが大きい。次にそれぞれのクラスタの平均値と分散は、
 $mean
          [,1]     [,2]     [,3]     [,4]     [,5]    [,6]     [,7]
 [1,] 672.5915 550.8010 350.8723 348.1141 721.9350 609.336 477.4751
 [2,] 652.4011 578.4403 238.2852 560.2122 701.1421 581.061 466.0265
             [,8]     [,9]
 [1,]   0.3521280 903.9811
 [2,] 252.0523710 941.5063
 
 $variance$sigma
 , , 1
          [,1]     [,2]
 [1,] 4510.760 2559.932
 [2,] 2559.932 1932.424
 , , 2
          [,1]     [,2]
 [1,] 5840.583 2278.554
 [2,] 2278.554 1350.294
 , , 3
             [,1]       [,2]
 [1,] 18864.03112  -47.12685
 [2,]   -47.12685 1101.31457
 , , 4
          [,1]     [,2]
 [1,] 120829.0 57238.70
 [2,]  57238.7 37994.28
 , , 5
          [,1]     [,2]
 [1,] 3609.556 2464.675
 [2,] 2464.675 2164.224
 , , 6
          [,1]     [,2]
 [1,] 2413.461 1707.167
 [2,] 1707.167 1546.101
 , , 7
          [,1]     [,2]
 [1,] 15796.04 14299.09
 [2,] 14299.09 16354.06
 , , 8
            [,1]         [,2]
 [1,] 84.9818781   -0.9859035
 [2,] -0.9859035 1455.7868607
 , , 9
          [,1]     [,2]
 [1,] 47.08964 19.93518
 [2,] 19.93518 12.34935
となっている。大きいクラスタとしてクラスタ2(551,578)、クラスタ6(609,581)、クラスタ5の(722,701)を見ると、分布図(x-y座標とも、数値を1/4にしている=1024値を256等分、また前回の処理の都合上FSCの100以下、SSCの50以下を切り取ってしまっていることに注意、クラスタ推定は切り取り前のデータを利用)上の点(138,145)、 (152,145)、 (181,175)のピーク(前者は赤い点2つ、最後の点は右上に伸びた塊の頂上)と対応付けられる。

// &ref(graph_FSC-SSC_100113_24h_40C_+20-1_limit_csv_256_800x800.png);

この解析では、それぞれの分布に属する点の数を、モデルの形から推定できていることであり、具体的には上記のclassification tableに相当する。6000点の中からクラスタ5(右上に伸びているクラスタ)に1013点が属していることが分かる。

その他の出力は
 $variance$scale
 [1]  1470.86706  1641.55235  4557.74192 36256.59046  1318.05373   903.89915
 [7]  7339.30971   351.73076    13.56889
 
 $variance$shape
 [1] 4.1389259 0.2416086
 
 $variance$orientation
 , , 1
           [,1]       [,2]
 [1,] 0.8514047 -0.5245093
 [2,] 0.5245093  0.8514047
 , , 2
           [,1]       [,2]
 [1,] 0.9224601  0.3860924
 [2,] 0.3860924 -0.9224601
 , , 3
              [,1]         [,2]
 [1,] -0.999996481 -0.002653105
 [2,]  0.002653105 -0.999996481
 , , 4
            [,1]       [,2]
 [1,] -0.8905668 -0.4548525
 [2,] -0.4548525  0.8905668
 , , 5
           [,1]       [,2]
 [1,] 0.8004262  0.5994314
 [2,] 0.5994314 -0.8004262
 , , 6
           [,1]       [,2]
 [1,] 0.7893714  0.6139159
 [2,] 0.6139159 -0.7893714
 , , 7
           [,1]       [,2]
 [1,] 0.7001755  0.7139708
 [2,] 0.7139708 -0.7001755
 , , 8
               [,1]          [,2]
 [1,] -0.0007192144 -0.9999997414
 [2,]  0.9999997414 -0.0007192144
 , , 9
            [,1]       [,2]
 [1,] -0.9102023 -0.4141639
 [2,] -0.4141639  0.9102023

トップ   編集 差分 バックアップ 添付 複製 名前変更 リロード   新規 一覧 単語検索 最終更新   ヘルプ   最終更新のRSS