Lecture10: 類似度計算, 階層的クラスター分析

前回の宿題+補足: app_lec09_ans

# create a sliter with the specified range
    sliderInput(inputId ="rank", 
                label ="X axis Range",
                min = 1, max =100, value=c(1,20), step=1),...

#Cf. app_pch1            
    sliderInput(inputId ="setCex", 
              label ="Marker Size (cex: character expansion)",
              min = 1, max =15, value=2, step=2)

類似度計算

  • 準備
  • proxyパッケージのインストール
  • パッケージの読み込み
  • univ頻度表の作成

proxyパッケージのインストール

install.packages("proxy", dependencies = TRUE)

proxyパッケージの読み込み

library(proxy)

ディレクトリ” univ”から出現単語行列を作成

getFreqDir関数の読み込み

source("utils3.R")

univディレクトリ内の頻度表の作成

univTable <- getFreqDir("univ")
head(univTable)

相関係数

ピアソン積率相関係数

\[Corr(x,y)= \frac{\sum (x_{i}-\overline{x}) (y_{i}-\overline{y})}{\sqrt{\sum (x_{i}-\overline{x})^2\sum (y_{i}-\overline{y})^2}} \]

相関係数行列(テキスト間)

tf <- getFreqDir("univ")
res <-cor(tf)
round(res,2)
        kyoto1 kyoto2 osaka1 osaka2 osaka3 osaka4 tokyo1 tokyo2 waseda1 waseda2
kyoto1    1.00   0.85   0.76   0.83   0.87   0.84   0.81   0.82    0.86    0.81
kyoto2    0.85   1.00   0.76   0.83   0.85   0.81   0.83   0.78    0.82    0.75
osaka1    0.76   0.76   1.00   0.84   0.80   0.79   0.71   0.73    0.76    0.68
osaka2    0.83   0.83   0.84   1.00   0.89   0.84   0.80   0.78    0.80    0.72
osaka3    0.87   0.85   0.80   0.89   1.00   0.88   0.84   0.80    0.81    0.75
osaka4    0.84   0.81   0.79   0.84   0.88   1.00   0.81   0.82    0.81    0.76
tokyo1    0.81   0.83   0.71   0.80   0.84   0.81   1.00   0.84    0.77    0.73
tokyo2    0.82   0.78   0.73   0.78   0.80   0.82   0.84   1.00    0.83    0.78
waseda1   0.86   0.82   0.76   0.80   0.81   0.81   0.77   0.83    1.00    0.86
waseda2   0.81   0.75   0.68   0.72   0.75   0.76   0.73   0.78    0.86    1.00

変数間の相関行列

転置(transpose): t関数

testMtx <- matrix(1:6, nrow=3)
testMtx
     [,1] [,2]
[1,]    1    4
[2,]    2    5
[3,]    3    6
# transpose
t(testMtx)
     [,1] [,2] [,3]
[1,]    1    2    3
[2,]    4    5    6

行と列を転置(transpose)する

corr <- simil(t(tf))
round(corr, 2)
        kyoto1 kyoto2 osaka1 osaka2 osaka3 osaka4 tokyo1 tokyo2 waseda1
kyoto2    0.85                                                         
osaka1    0.76   0.76                                                  
osaka2    0.83   0.83   0.84                                           
osaka3    0.87   0.85   0.80   0.89                                    
osaka4    0.84   0.81   0.79   0.84   0.88                             
tokyo1    0.81   0.83   0.71   0.80   0.84   0.81                      
tokyo2    0.82   0.78   0.73   0.78   0.80   0.82   0.84               
waseda1   0.86   0.82   0.76   0.80   0.81   0.81   0.77   0.83        
waseda2   0.81   0.75   0.68   0.72   0.75   0.76   0.73   0.78    0.86

変数間の相関行列(対角行列(diagonal)で結果表示)

res.corr <- simil(t(tf), diag=T)
round(res.corr, 2)
        kyoto1 kyoto2 osaka1 osaka2 osaka3 osaka4 tokyo1 tokyo2 waseda1 waseda2
kyoto1      NA                                                                 
kyoto2    0.85     NA                                                          
osaka1    0.76   0.76     NA                                                   
osaka2    0.83   0.83   0.84     NA                                            
osaka3    0.87   0.85   0.80   0.89     NA                                     
osaka4    0.84   0.81   0.79   0.84   0.88     NA                              
tokyo1    0.81   0.83   0.71   0.80   0.84   0.81     NA                       
tokyo2    0.82   0.78   0.73   0.78   0.80   0.82   0.84     NA                
waseda1   0.86   0.82   0.76   0.80   0.81   0.81   0.77   0.83      NA        
waseda2   0.81   0.75   0.68   0.72   0.75   0.76   0.73   0.78    0.86      NA

テキスト間のコサイン類似度

\[Cos(x,y)= \frac{\sum x_{i} y_{i}}{\sqrt{\sum x_{i}^2\sum y_{i}^2}} \]

res.cos <-simil(t(tf), method="cosine")

結果の形式:Distant Object出力

round(res.cos,2)
        kyoto1 kyoto2 osaka1 osaka2 osaka3 osaka4 tokyo1 tokyo2 waseda1
kyoto2    0.86                                                         
osaka1    0.77   0.77                                                  
osaka2    0.83   0.84   0.84                                           
osaka3    0.88   0.85   0.81   0.90                                    
osaka4    0.85   0.82   0.80   0.85   0.88                             
tokyo1    0.81   0.83   0.72   0.80   0.84   0.82                      
tokyo2    0.83   0.79   0.74   0.79   0.80   0.83   0.84               
waseda1   0.87   0.83   0.77   0.81   0.82   0.82   0.77   0.84        
waseda2   0.82   0.76   0.70   0.73   0.76   0.77   0.73   0.80    0.87

結果の形式:Matrix出力

res2.cos<-simil(t(tf), method="cosine", diag=T)
res2.cos<-as.matrix(res2.cos)
res2.cos[is.na(res2.cos)] <- 1
round(res2.cos,2)
        kyoto1 kyoto2 osaka1 osaka2 osaka3 osaka4 tokyo1 tokyo2 waseda1 waseda2
kyoto1    1.00   0.86   0.77   0.83   0.88   0.85   0.81   0.83    0.87    0.82
kyoto2    0.86   1.00   0.77   0.84   0.85   0.82   0.83   0.79    0.83    0.76
osaka1    0.77   0.77   1.00   0.84   0.81   0.80   0.72   0.74    0.77    0.70
osaka2    0.83   0.84   0.84   1.00   0.90   0.85   0.80   0.79    0.81    0.73
osaka3    0.88   0.85   0.81   0.90   1.00   0.88   0.84   0.80    0.82    0.76
osaka4    0.85   0.82   0.80   0.85   0.88   1.00   0.82   0.83    0.82    0.77
tokyo1    0.81   0.83   0.72   0.80   0.84   0.82   1.00   0.84    0.77    0.73
tokyo2    0.83   0.79   0.74   0.79   0.80   0.83   0.84   1.00    0.84    0.80
waseda1   0.87   0.83   0.77   0.81   0.82   0.82   0.77   0.84    1.00    0.87
waseda2   0.82   0.76   0.70   0.73   0.76   0.77   0.73   0.80    0.87    1.00

階層的クラスター分析

複数テキストでの単語の出現頻度の比較をする場合

相対頻度行列

relative.univTable <- getFreqDir("univ", relative=TRUE)
head(relative.univTable)

階層的クラスター分析: デフォルト値:complete法, euclidean距離))

hc <- hclust(dist(t(relative.univTable)))
plot(hc)

#rect.hclust(hc, k=3, border="red")

階層的クラスター分析: ward法 & euclidean距離

hc <- hclust(dist(t(relative.univTable)), method = "ward.D2")
plot(hc)

#rect.hclust(hc, k=3, border="red")

階層的クラスター分析: ward法 & キャンベラ距離

hc <- hclust(dist(t(relative.univTable), method = "canberra"), method = "ward.D2")
plot(hc)

#rect.hclust(hc, k=3, border="red")

課題3(締め切り2023年1月11日)

“app_asgnmt3”を編集し、Zipfの理論計算値の定数AとK、ならびに、頻度数プロットマーカーのタイプをインタラクティブに変更できるように拡張しなさい。

作成例: Shiny apps Demo

UIの情報

  • 定数K: 最小値=10, 最大値=100, 初期値=50
  • 定数A: 最小値=0.5, 最大値=1.5, 初期値=0.8, ステップ=0.1
  • プロットマーカーのタイプ: “app_pch2”の例を参考
LS0tCnRpdGxlOiAiTGVjMTA6IChGYWxsIDIwMjIpIgpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sKLS0tCgojIExlY3R1cmUxMDog6aGe5Ly85bqm6KiI566XLCDpmo7lsaTnmoTjgq/jg6njgrnjgr/jg7zliIbmnpAKCiMjIOWJjeWbnuOBruWuv+mhjCvoo5zotrM6IGFwcF9sZWMwOV9hbnMKKiA8YSBocmVmPSJodHRwczovL3JzdHVkaW8uZ2l0aHViLmlvL0RULyIgdGFyZ2V0PSJfYmxhbmsiPkRUIHBhY2thZ2U8L2E+CiogPGEgaHJlZj0iaHR0cHM6Ly9zaGlueS5yc3R1ZGlvLmNvbS9hcnRpY2xlcy9jc3MuaHRtbCIgdGFyZ2V0PSJfYmxhbmsiPnd3dyBzdWJmb2xkZXI8L2E+CiogdWkKYGBgCiMgY3JlYXRlIGEgc2xpdGVyIHdpdGggdGhlIHNwZWNpZmllZCByYW5nZQogICAgc2xpZGVySW5wdXQoaW5wdXRJZCA9InJhbmsiLCAKICAgICAgICAgICAgICAgIGxhYmVsID0iWCBheGlzIFJhbmdlIiwKICAgICAgICAgICAgICAgIG1pbiA9IDEsIG1heCA9MTAwLCB2YWx1ZT1jKDEsMjApLCBzdGVwPTEpLC4uLgoKI0NmLiBhcHBfcGNoMSAgICAgICAgICAgIAogICAgc2xpZGVySW5wdXQoaW5wdXRJZCA9InNldENleCIsIAogICAgICAgICAgICAgIGxhYmVsID0iTWFya2VyIFNpemUgKGNleDogY2hhcmFjdGVyIGV4cGFuc2lvbikiLAogICAgICAgICAgICAgIG1pbiA9IDEsIG1heCA9MTUsIHZhbHVlPTIsIHN0ZXA9MikKYGBgCgojIyDpoZ7kvLzluqboqIjnrpcKKiDmupblgpkKKyA8YSBocmVmPSJodHRwczovL3VyaWJvLmdpdGh1Yi5pby9ycGtnX3Nob3djYXNlL2RhdGEtYW5hbHlzaXMvcHJveHkuaHRtbCIgdGFyZ2V0PSJfYmxhbmsiPnByb3h5PC9hPuODkeODg+OCseODvOOCuOOBruOCpOODs+OCueODiOODvOODqworIOODkeODg+OCseODvOOCuOOBruiqreOBv+i+vOOBvworIHVuaXbpoLvluqbooajjga7kvZzmiJAKCiMjIyBwcm94eeODkeODg+OCseODvOOCuOOBruOCpOODs+OCueODiOODvOODqwpgYGB7ciwgZXZhbD1GQUxTRX0KaW5zdGFsbC5wYWNrYWdlcygicHJveHkiLCBkZXBlbmRlbmNpZXMgPSBUUlVFKQpgYGAKCiMjIyBwcm94eeODkeODg+OCseODvOOCuOOBruiqreOBv+i+vOOBvwpgYGB7ciwgZXZhbD1GQUxTRX0KbGlicmFyeShwcm94eSkKYGBgCgojIyDjg4fjgqPjg6zjgq/jg4jjg6oiIHVuaXYi44GL44KJ5Ye654++5Y2Y6Kqe6KGM5YiX44KS5L2c5oiQCiMjICBnZXRGcmVxRGly6Zai5pWw44Gu6Kqt44G/6L6844G/CmBgYHtyfQpzb3VyY2UoInV0aWxzMy5SIikKYGBgCgojIHVuaXbjg4fjgqPjg6zjgq/jg4jjg6rlhoXjga7poLvluqbooajjga7kvZzmiJAKYGBge3J9CnVuaXZUYWJsZSA8LSBnZXRGcmVxRGlyKCJ1bml2IikKaGVhZCh1bml2VGFibGUpCmBgYAoKIyMg55u46Zai5L+C5pWwCiMjIOODlOOCouOCveODs+epjeeOh+ebuOmWouS/guaVsAokJENvcnIoeCx5KT0gXGZyYWN7XHN1bSAoeF97aX0tXG92ZXJsaW5le3h9KSAoeV97aX0tXG92ZXJsaW5le3l9KX17XHNxcnR7XHN1bSAoeF97aX0tXG92ZXJsaW5le3h9KV4yXHN1bSAoeV97aX0tXG92ZXJsaW5le3l9KV4yfX0gJCQKCiMjIyDnm7jplqLkv4LmlbDooYzliJfvvIjjg4bjgq3jgrnjg4jplpPvvIkKYGBge3J9CnRmIDwtIGdldEZyZXFEaXIoInVuaXYiKQpyZXMgPC1jb3IodGYpCnJvdW5kKHJlcywyKQpgYGAKCiMjIyDlpInmlbDplpPjga7nm7jplqLooYzliJcKIyMjIyDou6Lnva7vvIh0cmFuc3Bvc2XvvIk6IHTplqLmlbAKYGBge3J9CnRlc3RNdHggPC0gbWF0cml4KDE6NiwgbnJvdz0zKQp0ZXN0TXR4CiMgdHJhbnNwb3NlCnQodGVzdE10eCkKYGBgCgojIyMjIOihjOOBqOWIl+OCkui7oue9ru+8iHRyYW5zcG9zZe+8ieOBmeOCiwpgYGB7cn0KY29yciA8LSBzaW1pbCh0KHRmKSkKcm91bmQoY29yciwgMikKYGBgCgojIyMg5aSJ5pWw6ZaT44Gu55u46Zai6KGM5YiX77yI5a++6KeS6KGM5YiXKGRpYWdvbmFsKeOBp+e1kOaenOihqOekuu+8iQpgYGB7cn0KcmVzLmNvcnIgPC0gc2ltaWwodCh0ZiksIGRpYWc9VCkKcm91bmQocmVzLmNvcnIsIDIpCmBgYAoKIyMg44OG44Kt44K544OI6ZaT44Gu44Kz44K144Kk44Oz6aGe5Ly85bqmCiQkQ29zKHgseSk9IFxmcmFje1xzdW0geF97aX0geV97aX19e1xzcXJ0e1xzdW0geF97aX1eMlxzdW0geV97aX1eMn19ICQkCgpgYGB7cn0KcmVzLmNvcyA8LXNpbWlsKHQodGYpLCBtZXRob2Q9ImNvc2luZSIpCmBgYAoKIyMjIOe1kOaenOOBruW9ouW8j++8mkRpc3RhbnQgT2JqZWN05Ye65YqbCmBgYHtyfQpyb3VuZChyZXMuY29zLDIpCmBgYAoKIyMjIOe1kOaenOOBruW9ouW8j++8mk1hdHJpeOWHuuWKmwpgYGB7cn0KcmVzMi5jb3M8LXNpbWlsKHQodGYpLCBtZXRob2Q9ImNvc2luZSIsIGRpYWc9VCkKcmVzMi5jb3M8LWFzLm1hdHJpeChyZXMyLmNvcykKcmVzMi5jb3NbaXMubmEocmVzMi5jb3MpXSA8LSAxCnJvdW5kKHJlczIuY29zLDIpCmBgYAoKIyMg6ZqO5bGk55qE44Kv44Op44K544K/44O85YiG5p6QCiogPGEgaHJlZj0iaHR0cHM6Ly93d3cuYWxiZXJ0MjAwNS5jby5qcC9rbm93bGVkZ2UvZGF0YV9taW5pbmcvY2x1c3Rlci9jbHVzdGVyX3N1bW1hcnkiIHRhcmdldD0iX2JsYW5rIj7lj4LogIPjgrXjgqTjg4gxPC9hPgoqIDxhIGhyZWY9Imh0dHBzOi8vd3d3MS5kb3NoaXNoYS5hYy5qcC9+bWppbi9SL0NoYXBfMjgvMjguaHRtbCIgdGFyZ2V0PSJfYmxhbmsiPuWPguiAg+OCteOCpOODiDI8L2E+CgojIyMg6KSH5pWw44OG44Kt44K544OI44Gn44Gu5Y2Y6Kqe44Gu5Ye654++6aC75bqm44Gu5q+U6LyD44KS44GZ44KL5aC05ZCICiMjIyMg55u45a++6aC75bqm6KGM5YiXCmBgYHtyfQpyZWxhdGl2ZS51bml2VGFibGUgPC0gZ2V0RnJlcURpcigidW5pdiIsIHJlbGF0aXZlPVRSVUUpCmhlYWQocmVsYXRpdmUudW5pdlRhYmxlKQpgYGAKCgojIyMg6ZqO5bGk55qE44Kv44Op44K544K/44O85YiG5p6QOiDjg4fjg5Xjgqnjg6vjg4jlgKTvvJpjb21wbGV0ZeazlSwgZXVjbGlkZWFu6Led6ZuiKSkKYGBge3J9CmhjIDwtIGhjbHVzdChkaXN0KHQocmVsYXRpdmUudW5pdlRhYmxlKSkpCnBsb3QoaGMpCiNyZWN0LmhjbHVzdChoYywgaz0zLCBib3JkZXI9InJlZCIpCmBgYAojIyMg6ZqO5bGk55qE44Kv44Op44K544K/44O85YiG5p6QOiB3YXJk5rOVICYgZXVjbGlkZWFu6Led6ZuiCmBgYHtyfQpoYyA8LSBoY2x1c3QoZGlzdCh0KHJlbGF0aXZlLnVuaXZUYWJsZSkpLCBtZXRob2QgPSAid2FyZC5EMiIpCnBsb3QoaGMpCiNyZWN0LmhjbHVzdChoYywgaz0zLCBib3JkZXI9InJlZCIpCmBgYAoKIyMjIOmajuWxpOeahOOCr+ODqeOCueOCv+ODvOWIhuaekDogd2FyZOazlSAmIOOCreODo+ODs+ODmeODqei3nembogpgYGB7cn0KaGMgPC0gaGNsdXN0KGRpc3QodChyZWxhdGl2ZS51bml2VGFibGUpLCBtZXRob2QgPSAiY2FuYmVycmEiKSwgbWV0aG9kID0gIndhcmQuRDIiKQpwbG90KGhjKQojcmVjdC5oY2x1c3QoaGMsIGs9MywgYm9yZGVyPSJyZWQiKQpgYGAKCiMjICDoqrLpoYwz77yI57eg44KB5YiH44KKMjAyM+W5tDHmnIgxMeaXpe+8iQojIyMgImFwcF9hc2dubXQzIuOCkue3qOmbhuOBl+OAgVppcGbjga7nkIboq5boqIjnrpflgKTjga7lrprmlbBB44GoS+OAgeOBquOCieOBs+OBq+OAgemgu+W6puaVsOODl+ODreODg+ODiOODnuODvOOCq+ODvOOBruOCv+OCpOODl+OCkuOCpOODs+OCv+ODqeOCr+ODhuOCo+ODluOBq+WkieabtOOBp+OBjeOCi+OCiOOBhuOBq+aLoeW8teOBl+OBquOBleOBhOOAggoKIyMjIOS9nOaIkOS+izogPGEgaHJlZj0iaHR0cHM6Ly9jb3B1bGFiby5zaGlueWFwcHMuaW8vYXBwX2FzZ210M19hbnMvIiB0YXJnZXQ9Il9ibGFuayI+U2hpbnkgYXBwcyBEZW1vPC9hPgojIyMjIFVJ44Gu5oOF5aCxIAoqIOWumuaVsEs6IOacgOWwj+WApD0xMCwg5pyA5aSn5YCkPTEwMCwg5Yid5pyf5YCkPTUwCiog5a6a5pWwQTog5pyA5bCP5YCkPTAuNSwg5pyA5aSn5YCkPTEuNSwg5Yid5pyf5YCkPTAuOCwg44K544OG44OD44OXPTAuMQoqIOODl+ODreODg+ODiOODnuODvOOCq+ODvOOBruOCv+OCpOODlzogImFwcF9wY2gyIuOBruS+i+OCkuWPguiAgwoK