思考ノイズ

無い知恵を絞りだす。無理はしない。

<はじパタ> 10.3 階層型クラスタリング (融合型)

K-平均法(10.2)を扱った前回に引き続きまして、クラスタリングの10章。今回は階層型、融合型というような表現がされている分類方法についてみていきたいと思います。

階層型クラスタリング

この融合型の方法は「定義」により距離的に近いとされる2つのクラスを一つのクラスとして融合することを繰り返していきます。今回は引き続きアヤメのデータを使いましたが、この場合はクラスは3つとなりますので、融合をすすめていき、最終的に3つのクラスに収束したら終わりとしました。

融合するための距離の4つの「定義」

融合にはクラス間の距離的な近さを計測します。はじパタでは距離の定義について以下の4つを紹介しています。

  • 単連結法

2クラスに属する要素通しで最も近いもの通しをそのクラス間の近さと定義して、その値が小さいもの通しを融合していく方法

  • 完全連結法

2クラスに属する要素通しで最も遠いもの通しをそのクラス間の近さと定義して、その値が小さいもの通しを融合していく方法

  • 群平均法

クラスタ内のすべての要素の平均をとり、各クラスの平均値間の距離をクラスの近さと定義する。その値が小さいもの通しを融合していく方法

  • ウォード法 クラスタを融合した場合の平均値が、もとのそれぞれのクラスタ内の平均値から変化する量をクラスタの近さと定義する。その値の小さいもの通しを融合していく方法

コーディングの際の考慮すべき点

アヤメデータは150の要素から成り立つのでそれぞれの要素の距離を示す行列を持つ必要があり、150 x 150のマトリックスとなります。さらにそのうち、上記の定義により近いと判断された2点が同じクラスタとされていきます。要素に番号を振り、融合された要素をリスト化していく形で、マトリックスが縮小(=クラスタ(リスト)内の要素が増えていく)ようにして、最終的に3つのクラスタになった時点で止めるようにしました。 最小値をもとめるコードの大部分は使いまわしが可能で、変更すべき点は、マトリックスのマージのプロセスのみとなりました。

実行結果

アヤメデータ4つの特性をそれぞれ2つずつ、計6通りのセットで確認しています。評価については10.2のK-平均法と同じく正解率とF値を出すようにしています。詳細はこちら

 ===  sepal length (cm)  ===  sepal width (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.35333333333333333
   * F-Value     --  0.502442878626
 ***  Complete Linkage Method
   * correctness --  0.52
   * F-Value     --  0.518563070642
 ***  Group Average Method
   * correctness --  0.7133333333333334
   * F-Value     --  0.741125961406
 ***  Ward Method
   * correctness --  0.78
   * F-Value     --  0.73182865776
 ===  sepal length (cm)  ===  petal length (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.6733333333333333
   * F-Value     --  0.763119756663
 ***  Complete Linkage Method
   * correctness --  0.84
   * F-Value     --  0.825008207203
 ***  Group Average Method
   * correctness --  0.7466666666666667
   * F-Value     --  0.758906000869
 ***  Ward Method
   * correctness --  0.84
   * F-Value     --  0.825008207203
 ===  sepal length (cm)  ===  petal width (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.34
   * F-Value     --  0.501625617698
 ***  Complete Linkage Method
   * correctness --  0.5866666666666667
   * F-Value     --  0.648999234975
 ***  Group Average Method
   * correctness --  0.7066666666666667
   * F-Value     --  0.745888252968
 ***  Ward Method
   * correctness --  0.8266666666666667
   * F-Value     --  0.806976194904
 ===  sepal width (cm)  ===  petal length (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.66
   * F-Value     --  0.760564623487
 ***  Complete Linkage Method
   * correctness --  0.94
   * F-Value     --  0.941265889017
 ***  Group Average Method
   * correctness --  0.8733333333333333
   * F-Value     --  0.863216837848
 ***  Ward Method
   * correctness --  0.8666666666666667
   * F-Value     --  0.854987973338
 ===  sepal width (cm)  ===  petal width (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.66
   * F-Value     --  0.760564623487
 ***  Complete Linkage Method
   * correctness --  0.8066666666666666
   * F-Value     --  0.738999798162
 ***  Group Average Method
   * correctness --  0.66
   * F-Value     --  0.760564623487
 ***  Ward Method
   * correctness --  0.84
   * F-Value     --  0.82020552181
 ===  petal length (cm)  ===  petal width (cm)  === 
 ***  Single Linkage Method
   * correctness --  0.6733333333333333
   * F-Value     --  0.763119756663
 ***  Complete Linkage Method
   * correctness --  0.96
   * F-Value     --  0.960162654599
 ***  Group Average Method
   * correctness --  0.96
   * F-Value     --  0.960162654599
 ***  Ward Method
   * correctness --  0.8933333333333333
   * F-Value     --  0.887000529608

考察

全体的な傾向としては単連結法が比較して正解率が良くなく、ほかの方法はケースバイケースでしょうか。そのなかでもウォード法は安定して高い値がでているように見えます。

K-平均法の時と同じようにすべての特性を扱った4次元データでの精度を確認してみましたが、やはり特徴点が顕著となる2点を使ったデータのほうが精度が高いという結果になりました。

また、K-平均法のほうが精度は高いです。結構苦労してコーディングしたのに悲しい限り。アヤメデータ以外でのK-平均法よりよくなるケースがあるのか気になるところではあります。

ソースコード

いかJupyter notebook上で動くソースです。

from sklearn import datasets
import matplotlib.pyplot as plt
from collections import Counter
#inline pyplot

class ClastaringFusion():
    def __init__(self):
        self.x = []
        self.y = []
        self.p = []
        self.t = []
        self.pt = []
        self.num = [50, 50, 50]
        iris = datasets.load_iris()
        
        self.feature_names = iris.feature_names
        features = iris.data
        targets = iris.target
        self.dmatrix = np.array([])

    def install_val(self,a=1,b=2):
        self.p = [a, b]
        self.x = features.T[self.p[0]]
        self.y = features.T[self.p[1]]
        self.t = targets
        self.distance()

        
    def distance(self, x=None, y=None):
        if x==None: x = self.x

        self.dmatrix = np.array([[np.sqrt((self.x[i] - self.x[j])**2 + (self.y[i] - self.y[j])**2 ) for i in range(150)] for j in range(150)]) 
        
        #return np.array([[np.sqrt((i - j)**2 + (i - j)**2) for i in x] for j in y])          
    def comp_fusion(self, list1, list2, fm=0):
        def ave(lst):
            x_lst = [self.x[i] for i in lst]
            y_lst = [self.y[i] for i in lst]

            m_x = sum(x_lst)/len(x_lst)
            m_y = sum(y_lst)/len(y_lst)


            return sum([(x_lst[i] - m_x)**2 + (y_lst[i] - m_y)**2 for i in range(len(lst))])
            
        if fm ==0:
            v = min([min([self.dmatrix[l1, l2] for l2 in list2]) for l1 in list1])
        elif fm==1:
            v =  max([max([self.dmatrix[l1, l2] for l2 in list2]) for l1 in list1])
        elif fm==2:
            v = sum([sum([self.dmatrix[l1, l2] for l2 in list2]) for l1 in list1]) /(len(list1) * len(list2))
        elif fm==3:
            v = ave(list1+list2) - (ave(list1) + ave(list2))    
        else:
            print("Error: Invalid number for Fusion Method")
            sys.exit(1)
            
        return v
    
    def fusion(self, fm=0, tree=None):
        '''
        fm (fusion_method)
        0: Single Linkage Method
        1: Complete Linkage Method
        2: Group Average Method
        3: Ward Method
        
        '''
        if tree==None:
            tree = [[i] for i in range(150)]
        dmtx = np.array([[self.comp_fusion(i, j, fm) for i in tree] for j in tree])

        for i in range(dmtx.shape[0]):
            dmtx[i, i] = 100
        
        val = np.argmin(dmtx)
        #print(val)
        #print(dmtx.shape)
        #print("D", dmtx[val//dmtx.shape[0], val%dmtx.shape[0]])
        minval = np.array(np.where(dmtx==dmtx[val//dmtx.shape[0], val%dmtx.shape[0]]))

        ctree = []
        checked = []
        lgth = minval.shape[1]
        #print(minval)
        #print("lgth: ", lgth)
        for i in range(lgth): 
            flag = True
            one = tree[minval[0][i]]
            two = tree[minval[1][i]]
            checked += (one + two)
            inc_num = []
            
            #print("one+two", one, two)
            #print(" ctree ", ctree)
            if minval[0][i] != minval[1][i]:

                for j, t in enumerate(ctree):
                    if set(t) & set(one) or set(t) & set(two):
                        inc_num.append(j)
                #print(" ", inc_num, " ==== ", one, two)
                #print(" ", ctree)
                if len(inc_num) == 0:
                    ctree.append(one+two)
                else:
                    ctree[inc_num[0]] = list(set(one+two+ctree[inc_num[0]]))
                    if len(inc_num) > 1:
                        for n in inc_num:
                            ctree[inc_num[0]] = list(set(ctree[inc_num[0]]+ctree[n]))
                            
                        for k, n in enumerate(inc_num[::-1]):
                            if k == len(inc_num)-1:
                                break
                            ctree.pop(n)
                    
          
        for i in tree:
            if not set(i) & set(checked):
                ctree.append(i)
        #print(ctree)
        test = []

        for t in ctree:
            test += t
        
        missing =  (set(range(150)) - set(test))
        tomatch =  (set(test) - set(range(150)))
        #print("Missing", missing)
        #print("Too Match", tomatch)
        return ctree  
    def show_all_result(self):
        fm_string = [
            "Single Linkage Method",
            "Complete Linkage Method",
            "Group Average Method",
            "Ward Method"
        ]
        for i, (a, b) in enumerate(itertools.combinations(range(4), 2)):
            
            print(" === ", self.feature_names[a], " === ", self.feature_names[b], " === ")
            self.install_val(a, b)
            #self.plot_map()

            for fm in range(4):
                tree = cf.fusion(fm)
                c = 0
                while len(tree) > 3:
                    #print("======")
                    tree = cf.fusion(fm, tree)
                    #print(len(tree))
                    #if len(tree) < 4:
                    #    break

                    c += 1
                crct, fv = cf.correct_rate(tree, False)
                print(" *** ", fm_string[fm])
                print("   * correctness -- ", crct)
                print("   * F-Value     -- ", fv)
             
        
    def correct_rate(self, tree, verbose=True):

        trees = []
        crcts = 0.0
        fvs = 0.0
        

        for p in list(itertools.permutations(range(3))):

            pt_tree = []
            for i in range(150):
                for j in range(3):
                    if i in tree[j]:
                        pt_tree.append(p.index(j))
                        break
                        
            trees.append(pt_tree)
            #print("tree", pt_tree)
                
            cnt = np.array([[Counter(pt_tree[i*50:i*50+50])[j] for i in range(3)] for j in range(3)])
            odr = [np.argmax([Counter(pt_tree[j*50:j*50+50])[i] for i in range(3)]) for j in range(3)]
            crct = []
            crctnum = 0

            for i in range(150):
                if i // 50 == pt_tree[i]:
                    crctnum += 1

            fv = 1.0

            for i in range(3):

                x = np.count_nonzero(pt_tree[i*50:i*50+50] == odr[i])
                z = np.count_nonzero(pt_tree[(i+1)%3*50:(i+1)%3*50+50] == odr[i]) + \
                    np.count_nonzero(pt_tree[(i+2)%3*50:(i+2)%3*50+50] == odr[i])
                rec = x/50
                pre = 50/(50+z)
                fv *= 2*rec*pre/(rec+pre)

            if crctnum > crcts:
                fvs = fv
                crcts = crctnum
            #print("crct", crctnum)
        if verbose:
            print("CORRECTNESS: ", crcts/150)
            print("F-Value    : ", np.cbrt(fvs))
        return crcts/150, np.cbrt(fvs)
        
cf = ClastaringFusion()
cf.install_val(1, 2)
cf.distance()
cf.show_all_result()

関連記事

bython-chogo.hatenablog.com

bython-chogo.hatenablog.com

<はじパタ> 10.2 K-平均法 おまけ

今回は以下の記事のおまけになります。

bython-chogo.hatenablog.com

K-平均法はランダムにK個クラスタした平均値をとって、その平均値に近い距離にいる要素をそのクラスタに分類する、という試行を収束するまで繰り返すものでした。

じゃあ実際にその平均値は各試行ごとにどのように変化するのかというのをアニメーションでプロットしてみました。 それが以下です。四角いポイントが各色の平均値となり、繰り返しによるクラスタの変化を各ステップごとにアニメーションで示しています。

f:id:bython-chogo:20180403235608g:plain

なるほど、最初にランダムにアサインされたデータから各テリトリーごとに収束していく様子が見やすくなりました。これは面白いですね。

ただ、ランダムでアクセスするので、以下のようにちょっとおかしな収束してしまうこともあるようです。特殊なパターンの収束をさけるために、試行は一回ではなく複数回やるのがよさそうです。

f:id:bython-chogo:20180403235600g:plain

さて、次回は10.3融合法です。この方法もアヤメデータを使った分類をしてみて、その正確さを評価してみたいと思います。

<はじパタ> 10.2 非階層型クラスタリング (K-平均法)

いろいろごたごたしている中でも少しづつはじパタすすめておりまして、とりあえず10章のストックができたのでこちらの備忘録を掲載します。

10章はクラスタリングで教師データなしの状態でパターンニングをおこなう手法になります。定番のアヤメデータをこのK-平均法でクラスタリングをおこなうとどうなるかをスクリプトを書いてみてみました。

アヤメデータ

このアヤメデータは3種50セット、合計150個のアヤメそれぞれの特性(がく片の幅、がく片の長さ、弁の幅および 弁の長さ)を記載しています。

例えば、がく片の幅と弁の長さでデータをプロットした場合の2次元データは以下のプロットで表現できます。

f:id:bython-chogo:20180403003147p:plain

つまり、今回の問いはこの150個のデータを教師データなしでとりあえず3種区別してみる、というものになります。 アヤメデータのプロットのスクリプトは以下のサイトを参考にさせていただいています。

Python: Iris データセットを Matplotlib で可視化してみる | CUBE SUGAR STORAGE

今回のスクリプトについて

今回は視覚化しやすいように4つの特性から2つをとったデータについてK平均法を使って3つのクラスタに分けてみます。4つのデータから2つとるので、6種のデータとなります。さらに4つのデータ全部を使ってどう変化するかもみてみたいと思います。

評価方法

評価には2つの指標を使います。ひとつは正解率、もうひとつは精度と再現率を評価できるF値を求めます。以下のリンクが詳しいかと思います。

抑えておきたい評価指標「正解率」「精度」「再現率」 - HELLO CYBERNETICS

K平均法

K-平均法は最初にそれぞれの値のクラスタをランダムに決定します。今回は3個のアヤメデータなので3つのクラスタにわけます。そしてその3つのクラスタの平均値のポイントを求め、またそれぞれの値で一番距離が近い平均値のクラスタに変更し直します。これを繰り返してデータが更新されなくなったら終わりとするものです。

実行結果

がく片の幅と弁の長さの結果

まずは先ほども例として挙げた「がく片の幅」と「弁の長さ」についてK-平均法の結果をプロットしてしてみます。

  • K-平均法の結果

f:id:bython-chogo:20180403005841p:plain

  • 実際のデータ

f:id:bython-chogo:20180403131050p:plain

  • 評価値
 ===  sepal width (cm) ,  petal length (cm)  === 
CORRECTNESS:  0.3867 ,   F-Value    :  0.4761
CORRECTNESS:  0.94 ,     F-Value    :  0.7417
CORRECTNESS:  0.8467 ,   F-Value    :  0.8337
CORRECTNESS:  0.9 ,  F-Value    :  0.897
CORRECTNESS:  0.9333 ,   F-Value    :  0.9328
CORRECTNESS:  0.9533 ,   F-Value    :  0.954
CORRECTNESS:  0.96 ,     F-Value    :  0.9607
CORRECTNESS:  0.94 ,     F-Value    :  0.9413
CORRECTNESS:  0.9267 ,   F-Value    :  0.9271

もとのデータが比較的わかりやすく分かれていることもあるので、正答率もわるくないと思われます。ランダムにはじめて平均値とその距離でこうもきれいにわかれるものなのかと少し驚きました。3つの平均値をみてみると最初は3つとも近い場所にいるのですが、繰り返していくうちに平均的なポジションに落ち着く様子がわかりました。これも今度アニメーションでプロットしてみたいと思います。

また正解率やF値はピークを越えて少しだけ下がるところで収束してしまうようです。実際は回答とすべき教師データはないので、収束まつしか判断のしようがなさそうです。

4つの特性から2つピックアップしてみる

では、4つの特性から2つピックアップした時のそれぞれのK-平均法の評価値はどのように変わるか見てみましょう。初期値がランダムなので評価ごとに値は変更しますが、とりあえず評価値の大小関係は変わらないと考えてよいです。

 ===  sepal length (cm) ,  sepal width (cm)  === 
CORRECTNESS:  0.4 ,  F-Value    :  0.4831
CORRECTNESS:  0.7867 ,   F-Value    :  0.689
CORRECTNESS:  0.8067 ,   F-Value    :  0.8134
CORRECTNESS:  0.8267 ,   F-Value    :  0.8363
CORRECTNESS:  0.82 ,     F-Value    :  0.8306

 ===  sepal length (cm) ,  petal length (cm)  === 
CORRECTNESS:  0.3867 ,   F-Value    :  0.4652
CORRECTNESS:  0.8267 ,   F-Value    :  0.7216
CORRECTNESS:  0.88 ,     F-Value    :  0.8794
CORRECTNESS:  0.9 ,  F-Value    :  0.9021
CORRECTNESS:  0.8867 ,   F-Value    :  0.8918
CORRECTNESS:  0.88 ,     F-Value    :  0.885
CORRECTNESS:  0.8867 ,   F-Value    :  0.8897
CORRECTNESS:  0.88 ,     F-Value    :  0.8824

 ===  sepal length (cm) ,  petal width (cm)  === 
CORRECTNESS:  0.3667 ,   F-Value    :  0.4585
CORRECTNESS:  0.8333 ,   F-Value    :  0.6973
CORRECTNESS:  0.8267 ,   F-Value    :  0.83
CORRECTNESS:  0.8267 ,   F-Value    :  0.8352
CORRECTNESS:  0.8333 ,   F-Value    :  0.8414

 ===  sepal width (cm) ,  petal length (cm)  === 
CORRECTNESS:  0.3867 ,   F-Value    :  0.4761
CORRECTNESS:  0.94 ,     F-Value    :  0.7417
CORRECTNESS:  0.8467 ,   F-Value    :  0.8337
CORRECTNESS:  0.9 ,  F-Value    :  0.897
CORRECTNESS:  0.9333 ,   F-Value    :  0.9328
CORRECTNESS:  0.9533 ,   F-Value    :  0.954
CORRECTNESS:  0.96 ,     F-Value    :  0.9607
CORRECTNESS:  0.94 ,     F-Value    :  0.9413
CORRECTNESS:  0.9267 ,   F-Value    :  0.9271

 ===  sepal width (cm) ,  petal width (cm)  === 
CORRECTNESS:  0.3933 ,   F-Value    :  0.4796
CORRECTNESS:  0.9533 ,   F-Value    :  0.7596
CORRECTNESS:  0.86 ,     F-Value    :  0.7473
CORRECTNESS:  0.8467 ,   F-Value    :  0.841
CORRECTNESS:  0.8733 ,   F-Value    :  0.876
CORRECTNESS:  0.8933 ,   F-Value    :  0.8965
CORRECTNESS:  0.9133 ,   F-Value    :  0.9168
CORRECTNESS:  0.9267 ,   F-Value    :  0.9293
CORRECTNESS:  0.9267 ,   F-Value    :  0.9292

 ===  petal length (cm) ,  petal width (cm)  === 
CORRECTNESS:  0.36 ,     F-Value    :  0.4523
CORRECTNESS:  0.9467 ,   F-Value    :  0.7468
CORRECTNESS:  0.8467 ,   F-Value    :  0.8292
CORRECTNESS:  0.9267 ,   F-Value    :  0.9258
CORRECTNESS:  0.96 ,     F-Value    :  0.9602
CORRECTNESS:  0.9733 ,   F-Value    :  0.9735
CORRECTNESS:  0.9667 ,   F-Value    :  0.9672
CORRECTNESS:  0.96 ,     F-Value    :  0.9606

概ね、0.85~0.95といったあたりに収まるでしょうか。

4つの特性を使った評価

では、2つといわず一気に4つの特性を使った4次元での評価はどうなるか確認をしてみたいと思います。

CORRECTNESS:  0.4267 ,    F-Value    :  0.5037
CORRECTNESS:  0.8467 ,   F-Value    :  0.7439
CORRECTNESS:  0.8667 ,   F-Value    :  0.8586
CORRECTNESS:  0.9 ,  F-Value    :  0.897
CORRECTNESS:  0.9 ,  F-Value    :  0.9007
CORRECTNESS:  0.92 ,     F-Value    :  0.9221
CORRECTNESS:  0.92 ,     F-Value    :  0.9227
CORRECTNESS:  0.9067 ,   F-Value    :  0.9096
CORRECTNESS:  0.9 ,  F-Value    :  0.9021
CORRECTNESS:  0.8933 ,   F-Value    :  0.8937
CORRECTNESS:  0.8867 ,   F-Value    :  0.8866

あれれ、精度としては特徴がよく表れている2つをピックアップしたほうがよい数値がでているようにもみえます。でもおそらく、その結果はその特徴量に強く適合がしすぎてしまっているかと思われるので、こちらの4次元での方法はより「一般的」「総合的」な評価がしやすくなるのかとは思います。

ソースコード

今回利用したソースコードです。Jupyter notebook上で動作を確認しており、 processで実行、2つの引数はアヤメの4つ特性 0~3のうち2つピックするのと、3つ目の引数はプロットをするかTrue/Falseで表しています。

import numpy as np
from sklearn import datasets
import matplotlib.pyplot as plt
from collections import Counter
#inline pyplot

class Iris_PM:
    def __init__(self):
        self.p = [2, 3]
        self.x = []
        self.y = []
        self.t = []
        self.pt = []
        self.num = [50, 50, 50]
        iris = datasets.load_iris()
        features = iris.data
        targets = iris.target
    def install(self, a=2, b=3):
        self.x = features.T[a]
        self.y = features.T[b]
        self.t = targets

    def plot_map(self, df=False):
        p = self.pt
        if df:
            p = self.t

        for t, marker, c in zip(range(3), '>ox', 'rgb'):
             
            plt.scatter(
                #self.x[50*t:50*t+49],
                #self.y[50*t:50*t+49],
                features[p == t, self.p[0]],
                features[p == t, self.p[1]],
                marker=marker,
                c=c,
            )
            plt.xlim(1.9, 4.5)
            plt.ylim(0.9, 7.1)
            plt.xlabel(feature_names[self.p[0]])
            plt.ylabel(feature_names[self.p[1]])
            plt.autoscale()
            plt.grid()
        #print(linear_pro(w, np.array([1.,1.,1.])))
        plt.show()
    def randomize(self):
        self.pt = np.array([np.random.randint(3) for i in range(sum(self.num))])
    
    def means(self):
        nu = [[0.0, 0.0], [0.0, 0.0], [0.0, 0.0]]
        cnt = [0, 0, 0]
        for i in range(sum(self.num)):
            nu[krs.pt[i]][0] += self.x[i]
            nu[krs.pt[i]][1] += self.y[i]
            cnt[krs.pt[i]] += 1
        for i in range(3):
            nu[i][0] /= cnt[i]
            nu[i][1] /= cnt[i]

        return nu
    def correct_rate(self):
        #print(self.pt)
        cnt = np.array([[Counter(self.pt[i*50:i*50+50])[j] for i in range(3)] for j in range(3)])
        #print(cnt)
        #print([np.argmax([Counter(self.pt[j*50:j*50+50])[i] for i in range(3)]) for j in range(3)])
        #print([np.argmax([Counter(self.pt[i*50:i*50+50])[j] for i in range(3)]) for j in range(3)])
        odr = [np.argmax([Counter(self.pt[j*50:j*50+50])[i] for i in range(3)]) for j in range(3)]
        lodr = []
        for i in odr:
            lodr += ([i]*50)

        fv = 1.0
        
        for i in range(3):
      
            x = np.count_nonzero(self.pt[i*50:i*50+50] == odr[i])
            z = np.count_nonzero(self.pt[(i+1)%3*50:(i+1)%3*50+50] == odr[i]) + \
                np.count_nonzero(self.pt[(i+2)%3*50:(i+2)%3*50+50] == odr[i])
            rec = x/50
            pre = 50/(50+z)
            fv *= 2*rec*pre/(rec+pre)

        print("CORRECTNESS: ", "{:.4}".format(np.sum(self.pt == lodr)/150), ",\t F-Value    : ", "{:.4}".format(np.cbrt(fv)))
        
    def process(self, a=1, b=2, plot_map=False):
        self.install(a,b)
        self.randomize()

        bak_list = []
        cnt = 0
        print(" === ", feature_names[a], ", ", feature_names[b]," === ")
        while self.pt is not bak_list and cnt < 30:

            self.correct_rate()

            cnt += 1
            bak_list = self.pt
            nu = self.means()

            up = []
            #krs.plot_map()
            for i in range(sum(self.num)):
                #print(n,":", krs.x[i], krs.y[i], nu)
                up.append(np.argmin([[(self.x[i]-n[0])**2 + (self.y[i]-n[1])**2 for n in nu]]))

            self.pt = np.array(up)
            nu = self.means()


            if np.allclose(bak_list, self.pt):
                if plot_map:
                    self.plot_map()
                    self.plot_map(True)
                break        
        print("")
krs = Iris_PM()

krs.process(1,2,False) 

CTS LabsのAMDプロセッサに対する脆弱性のPaperをまとめてみた

3月13日に公表されたAMD脆弱性について、イスラエルのセキュリティ会社、CTS Labsが出しているWhite Paperをまとめてみた。

AMDFLAWS

POC(Proof of Concept: コンセプトの 実効性検証)の開示もまだのようで具体的な手法の情報も乏しく、現時点で実現可能性に関しては疑問がのこるものになる。

いくつかのニュースサイトで記事になっているように、現時点での情報ではこの攻撃による深刻性については不明確だ。すでに強い権限を持った状態が前提の作業で、前段階のセキュリティが破られいることは条件になりそうだ。また、この手の脆弱性にはCert/CCから発行されるCVEナンバーがついて、公表されるのが定例だが、この会社はそのような公開プロセスを飛ばして公表に踏み切っているというため、お行儀がよくないと批判も受けているとのこと。 なんにせよ、最終的な結論は詳細情報のアップデートが必要そうだ。

CTSの以下の記述のとおり、現在、対応パッチや緩和策のためにAMD/Microsoftなどの関係企業でシェアしているとのこと。(3月15日現在)

Q. Doesn't this publication put users at risk?
A. No. All technical details that could be used to reproduce the vulnerabilities have been redacted from this publication. CTS has shared this information with AMD, Microsoft, and a small number of companies that could produce patches and mitigations.

以下、公式サイトおよび論文のまとめとなる。技術知識が深くない私が、Paperを読んで理解できた上っ面のまとめとなること、勘違い、知識不足をご容赦いただきたい。

脆弱性概要

CTS Labsによると今回公表されている脆弱性は大きく4つのクラス(MASTERKEY/RYZENFALL/FALLOUT/CHIMERA)に分けられ、13個の脆弱性があるという。

f:id:bython-chogo:20180314173350p:plain

公式サイトに記載されている各クラスの概況を大雑把にまとめる。

-- AMD Secure Processor firmwareに対する攻撃で、Firmware ベースのセキュリティ機能 (Secure Encrypted Virtualization, Firmware Trusted Platform Module)に対する書き換え、FirmwareベースでのRansomewareの書き込みが可能だという。

-- SMRAM(System Momory RAM), やWindowsが持っているセキュアなメモリ空間にアクセスすることが可能になるため、AMD Secure Processorをコントロールすることができる。MASTERKEY攻撃に必要な不正なBIOSコードの書き込みも行える。

-- SMRAM(System Momory RAM), やWindowsが持っているセキュアなメモリ空間にアクセスすることが可能になるため、AMD Secure Processorをコントロールすることができる。不正なBIOSコードの書き込みも行える。RYZENFALLはSecure OS機構に対する攻撃でClient/Workstation製品に対する攻撃だが、FALLOUTはBoot Loaderチップに対する攻撃で対象はEPYC Server製品となる

  • CHIMERA: Backdoors Inside Ryzen Chipset

-- Ryzenなどのプラットフォームで用いられている3rd Partyの隠れたBackdoorを突くもの。DMA(Direct Memory Access)エンジンなどに攻撃をすることで、プラットフォームのIOに対する不正攻撃が可能となる。Key Logger/Network Access/PypassMemory Protectionなどが挙げられている。

Paper まとめ

MASTERKEY: Unauthorized Code Execution and Malware Persistency on AMD Secure Processor

AMDのSecure Bootの機能にMalwareなど不正なコードをしこみ、最終的にSecure ProcessorのKernel Modeの挙動もコントロールが可能とされている。Malwareを仕込めば、Secure Bootを不正にBypassしてSecure Processorの機能を無効化もできる。 通常はSignatureでサインされているものでないとBIOSのアップデートはできないのだが、RYZENFALLもしくはFALLOUTでSystem Management Mode(SMM)という上級権限に入ることにより、それを可能としている。

RYZENFALL: Vulnerabilities in Ryzen Secure Processor

RyzenなどClient/Workstation向けのシステムに導入されているSecure OSはfTPM (Firmware Trusted Platform Moduel)に含まれる認証情報からのスタートとなる。このSecure BootのProcessで利用されるメモリ空間は特に権限が強い Fenced DRAMというところでSecure Processorでしか利用ができない空間になる。 ここに対するRYZENFALLの攻撃はadministrator 権限が必要となり、Secure Processorへのアクセスはベンダから公式に配布されるドライバを使用する。RYZENFALLに成功するとSecure Processorでの実行が認証なしで行えるようになり、ハードウェアで保護されているMemory Regionへのアクセスを可能とする。具体的には以下の攻撃チャネルがある。

  • Windows Isolated User Mode and Isolated Kernel Mode (VTL1)
  • Secure Management RAM (SMRAM)
  • AMD Secure Processor Fenced DRAM

これらを経由することで、 以下のことが実行可能になる。 - Microsoft Virtualization-based Securityをバイパスして Networkの権限を奪うことが可能。 - SMMにMalwareをいれることができOSやHypervisorでの攻撃が可能となる - BIOSのReflashに対するプロテクションの無効化 - VTL1にMalwareを入れることが可能。OS上でのエンドポイントでのSecurity機能の無効化が可能。 - AMD Secure Processor自身にMalwareの注入が可能 - fTPMなどのSecurity機能をバイパスしてのブートアップが可能

FALLOUT: Vulnerabilities in EPYC Server Secure Processor

FALLOUTはEPYCなどのサーバを対象とした攻撃方法。EPYCのSecurity Processorで利用するBoot loaderのコンポーネント内の脆弱性をつくと、Hardware Validated BootやSecure Encrypted VirtualizationなどのSecurityに影響を与えることができる。 RYZENFALLと同じく、administrator権限での実行が必要となり、公式のドライバを使用する。これを利用してHardwareで保護されているMemory Regionへのアクセスすることができる。具体的には以下の攻撃チャネルがある。

  • Windows Isolated User Mode and Isolated Kernel Mode (VTL1)
  • Secure Management RAM (SMRAM)

これらを経由することで、 以下のことが実行可能になる。 - Microsoft Virtualization-based Securityをバイパスして Networkの権限を奪うことが可能。 - SMMにMalwareをいれることができOSやHypervisorでの攻撃が可能となる - BIOSのReflashに対するプロテクションの無効化 - VTL1にMalwareを入れることが可能。OS上でのエンドポイントでのSecurity機能の無効化が可能。

CHIMERA: Backdoors Inside Ryzen Chipset

AMDのPromontory Chipset内にある、製造工程上などで利用されたとみられる隠れたバックドアを突く方法。2つのバックドアが存在している。一つはChipsetを実行するFirmwareにあり、ASICハードウェアハードの中にある。特に後者はChipset内にあるものなので修正が難しく、ワークアラウンドかRecallしかない。 Chipset内にあるので、提供されているSecurityの機能にも影響を及ぼすことが可能になり、さらにChipsetはUSB, SATA, PCI-EといったI/O周りも制御ができるので、LAN, Wifi, Bluetoothにも影響を及ぼすことができる。彼らの研究ではDMA(Direct Memory Access)に対する攻撃で、OSに影響を与えることができたという。このことにより、以下のような追加の攻撃ができるようになると考えられている。

  • Key Logger
  • Network Access
  • Bypass Memory Protection

<はじパタ> 6.1.1 超平面の方程式

はじめてのパターン認識 6.1.1 超平面の方程式に関してのメモです。

はじめてのパターン認識

はじめてのパターン認識

ここでは超境界面の数式  f(\mathbf{x}) = \mathbf{w} \mathbf{x} + w_0 を単位法線ベクトル  n と位置ベクトル P の関係に変形して、  f(\mathbf{x}) = \mathbf{n}^T (\mathbf{x} - \mathbf{P}) としている。

境界面へのベクトル P とそこからのびる、単位法線ベクトル  n を分離することにより、 n (x-P)内積の計算結果の正負を考えると、境界面のどちら側にいるかの識別が可能となる。 例題6.1の計算で  n = (\frac{2}{\sqrt{5}}, -\frac{1}{\sqrt{5}}), P = ( 1, 0) という値が計算された。

ここで、以下の3つの点について考えてみた。


a = (0, 3),
b = (2, 2),
c = (2, -1)

この三点と、ベクトル  n, P をプロットしてみる。

plt.axes().set_aspect('equal', 'datalim')

x = np.linspace(-1, 3)
y = 2*x - 2
plt.plot(x, y, "r-")

plt.plot(0, 3, "bo")
plt.text(0.2, 3, "a")

plt.plot(2, 2, "bo")
plt.text(2.2,2, "b")

plt.plot(2, -1, "bo")
plt.text(2.2,-1, "c")

plt.text(1.2,0, "P")
plt.text(2,-0.5, "n")

plt.quiver(0, 0, 1, 0, angles='xy',scale_units='xy',scale=1)
plt.quiver(1, 0, 2.0/np.sqrt(5), -1.0/np.sqrt(5), angles='xy',scale_units='xy',scale=1)

plt.grid()
plt.xlim(-2, 3)
plt.ylim(-2, 4)
plt.show()

f:id:bython-chogo:20180127235115p:plain

 \mathbf{n} (x - P)  x  a, b, c をそれぞれ代入すると以下のようになる。

 f(a) = \frac{2}{\sqrt{5}}  f(b) = 0  f(c) = \frac{-1}{\sqrt{5}}

この正負の結果は  (\mathbf{x}-\mathbf{P}) \mathbf{n} のつくる角度がそれぞれ、鋭角、直角、鈍角であるので、それぞれの結果は正、ゼロ、負となるわけである。

わかりやすくするために、 nベクトルを平行移動させて、図で表すと以下のようになる。

plt.axes().set_aspect('equal', 'datalim')

x = np.linspace(-1, 3)
y = 2*x - 2

plt.plot(0, 3, "bo")
plt.text(0.2, 3, "a")
plt.quiver(0, 3, 1, -3, angles='xy',scale_units='xy',scale=1)
plt.quiver(0, 3, 2.0/np.sqrt(5), -1.0/np.sqrt(5), angles='xy',scale_units='xy',scale=1)

plt.plot(2, 2, "bo")
plt.text(2.2,2, "b")
plt.quiver(2, 2, -1, -2, angles='xy',scale_units='xy',scale=1)
plt.quiver(2, 2, 2.0/np.sqrt(5), -1.0/np.sqrt(5), angles='xy',scale_units='xy',scale=1)

plt.plot(2, -1, "bo")
plt.text(2.2,-1, "c")
plt.quiver(2, -1, -1, 1, angles='xy',scale_units='xy',scale=1)
plt.quiver(2, -1, 2.0/np.sqrt(5), -1.0/np.sqrt(5), angles='xy',scale_units='xy',scale=1)

plt.text(1.2,0, "P")

plt.grid()
plt.xlim(-2, 3)
plt.ylim(-2, 4)
plt.show()

f:id:bython-chogo:20180127235323p:plain

お気づきかとは思いますが、正負が逆になってしまう。2クラスの識別なら正負が反転していても問題がないんだけど、気持ちが悪い。 本文の例題にもこっそり引き算を逆にして正負を反転させているような記述がみえるのだが、ごまかしているのか、私のベクトル計算の考え方がちがうのか判断がつかない。

機械学習のお勉強

久しぶりのエントリー。ドラフト版を書いては消去するのを繰り返しているうちに遠ざかってしまいました。

そもそも機械学習やプログラミングの勉強メモとして再開していたこのブログ、お休みの間も機械学習の勉強を続けていました。ただやっていくうちにどんどんとベーシックなほうへ向かっていき、最終的に機械学習(と統計)のための数学の本を買いだす始末。

ただこの本がとてもよかったです。統計、機械学習でつかう数学式のみを中1レベルから丁寧に教えて書いてくれていて、いかに自分が数学から離れて忘却されてしまっていたか、さらに機械学習本でわかったつもりになっていたかがよくわかりました。 またこれ以上Low Layerに向かうことはないだろう、ここからがスタートだと思うことができました。その後、一回、計算式の変換で挫折をした「はじパタ」を読み直しています。

はじめてのパターン認識

はじめてのパターン認識

Markdownで書いたメモもたまってきたので、ブログに放出をしていきたいと思います。

このほかの本も読んでいるのですが、それはまた機会があれば。

芸人のラジオ考察

またしても他人の褌で相撲を取るようなエントリーなのですが、思っちゃたんだからしょうがない。

hiko1985.hatenablog.com

芸人ラジオ大好きで主にニッポン放送TBSラジオの深夜枠はラジコを駆使しながら拾えるだけ拾うようにしていたのですが、ハライチのターンは存在を知りながらカバーしてませんでした。いわゆる食わず嫌いってやつです。しかしこの動画はそんな斜に構えた私の心を打ち砕くのに十二分なインパクトがありました。

脳内世界観強めの話芸

あぁ、ハライチのラジオは岩井さんのほうなんだなと。もちろん、澤部さんも面白いのですが、パンチがツボにクリーンヒットさせるのは岩井さんなんだと。岩井さんは独自の世界を脳内に構築するタイプの方で、それを現実世界に面白く脚色しながらアウトプットできる能力の持ち主なんでしょう。これは別名、心の闇ということもできますが。 たぶん、この手の脳内の世界っていうのは多かれ少なかれ誰しももっているもので、芸術関係の人はこの世界を音楽だったり絵画だったりに出力をするのですが、岩井さんは芸人としての話芸で話すタイプなんだと確信しました。類型でアルコアンドピースの平子さんもこうした世界観づくりは得意ですが、平子さんの場合はもっとファンタジーだったり中二病感の濃度が濃くなります。

芸人ラジオ勝手に世代分け

だんだんと自分の観測範囲内の自己研究が強めです。当然異論は認めます。

いまのTBS・ニッポン放送の芸人ラジオは世代で以下の3つに分けると考えやすいのかなと思っています。さらに代表的なラジオ番組も入れてみました。

20年戦士

10年戦士

5年戦士以内(前身番組なども含め)

新世代のラジオ芸人考察

ラジオは人気が大きくなると長寿になりやすく、安定したしゃべり手が長く番組を持つ傾向にあると思います。その中でしっかり新しい風を吹き込んでいって、それらがちょっと安定するかしないかのところが新鮮で面白い、そんな聴き時なところが5年戦士に入る番組なのかなと、思っていたりします。この人たちは同世代でもあるののでよくわかるのですが、アニメとか映画の表現が色濃くなった90年代に子供時代をすごし、かつ「大人になってもアニメを見てもいい」という新常識が生まれた境目にいた世代だとおものうので、最初に言ったような「脳内の世界観」の蓄積が充実してきている人が多いのかなと。 その脳内の世界観をデータベースとして、芸人としてのトーク力のスキルでうまく処理し、さらにラジオという特殊な場所という好条件が重なった結果、今回の岩井さんのような新しいラジオスターが生まれたのではと考えます。

今後も芸人ラジオから耳が離せませんです。

f:id:bython-chogo:20171104175934p:plain