「実践コンピュータビジョン」の全演習問題をやってみた。詳細編 第８章 画像認識

まえおき

先日投稿した「実践コンピュータビジョン」の全演習問題をやってみたの詳細です。なお、あくまで「やってみた」であって、模範解答ではありません。「やってみた」だけなので、うまくいかなかった、という結果もあります。

この記事は「実践コンピュータビジョン」（オライリージャパン、Jan Erik Solem著、相川愛三訳）の演習問題に関する記事です。この記事を書いた人間と書籍との間には何の関係もありません。著者のJan氏と訳者の相川氏にはコンピュータビジョンについて学ぶ機会を与えて頂いたことに感謝します。

実践コンピュータビジョンの原書("Programming Computer Vision with Python")の校正前の原稿は著者のJan Erik Solem氏がCreative Commons Licenseの元で公開されています。。

回答中で使用する画像の多く、また、sfm.py, sift.py, stereo.pyなどの教科書中で作られるファイルはこれらのページからダウンロードできます。
オライリージャパンのサポートページ
Jan Erik Solem氏のホームページ
相川愛三氏のホームページ

ここで扱うプログラムはすべてPython2.7上で動作を確認しています。一部の例外を除きほとんどがJupyter上で開発されました。また必要なライブラリーを適宜インストールしておく必要があります。Jupyter以外の点は、上記の書籍で説明されています。

結構な分量なので章毎に記事をアップしています。今回は第８章です。

８章　画像認識

最近では画像認識といえばディープラーニングな雰囲気ではありますが、トラディショナルなCVが役に立つ機会もまだまだ多いようです。近傍法、SVMなどによる画像認識を学びます。

8.5.1 The performance of the kNN classifier depends on the value of k. Try to vary this number and see how the accuracy changes. Plot the decision boundaries of the 2D point sets to see how they change.

k近傍法のkを変えて精度の変化を見ろ

私の回答

https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%201.ipynb
ひさしぶりに教科書の演習問題っぽくてほっとします。

from numpy import *
from PIL import *
import pickle
from pylab import *
import os
import sift
import dsift
import imtools

def read_gesture_features_labels(path):
    # make a list of the files with .dsift at the end
    featlist = [os.path.join(path, f) for f in os.listdir(path)
               if f.endswith('.dsift')]

    # read features
    features = []
    for featfile in featlist:
        l, d = sift.read_features_from_file(featfile)
        features.append(d.flatten())
    features = array(features)

    # generate labels
    labels = [featfile.split('/')[-1][0] for featfile in featlist]

    return features, array(labels)

features, labels = read_gesture_features_labels('train/')
test_features, test_labels = read_gesture_features_labels('test/')
classnames = unique(labels)

import knn
# k neighbors
accs = []
for k in arange(1, 10):
    knn_classifier = knn.KnnClassifier(labels, features)

    res = array([knn_classifier.classify(test_features[i], k) for i in
            range(len(test_labels))])

    # accuracy
    acc = sum(1.0*(res==test_labels))/len(test_labels)
    print 'Accuracy:', acc
    accs.append(acc)

結果
Accuracy: 0.716346153846
Accuracy: 0.605769230769
Accuracy: 0.519230769231
Accuracy: 0.475961538462
Accuracy: 0.427884615385
Accuracy: 0.365384615385
Accuracy: 0.331730769231
Accuracy: 0.288461538462
Accuracy: 0.269230769231

上から順にk=1, 2, ..., 10です。だんだん悪くなっているようです。

問題では「判別境界がどう変わるか観察しよう」というのもあります。回答はこちら
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%201-2.ipynb
結果の抜粋のみここにだしておきます。

k=1

k=2

k=3

中途略
k=10

8.5.2 The hand gesture data set in Figure 8.3 also contains images with more complex background (in the "complex/" folders). Try to train and test a classifier on these images. What is the difference in performance? Can you suggest improvements to the image descriptor?

図８−３の手話データと同梱された複雑な背景の画像で訓練とテストを行え

私の回答

やっていることは教科書の内容と同じなのでコードは省略します。結果は、

k=1
Accuracy: 0.26213592233

これはシンプルなデータセットでの結果とくらべるとだいぶ悪いですね。
シンプルな方の結果はこちらです。
k=1
Accuracy: 0.716346153846

8.5.3 Try to vary the number of dimensions after PCA projection of the gesture recognition features for the Bayes classifier. What is a good choice? Plot the singular values S, they should give a typical "knee" shaped curve as the one below. A good compromise between ability to generate the variability of the training data and keeping the number of dimensions low is usually found at a number before the

curve flattens out.
ベイズ分類器でPCA射影後に削減する次元数を変えろ

私の回答

これもやっていることは教科書の内容と大差ないのでコードは省きます。以下、手話データの複雑な方（train2）を使った結果です。kは削減後の次元数です。

k=12
Accuracy: 0.165048543689
k=25
Accuracy: 0.223300970874
k=50
Accuracy: 0.26213592233
k=75
Accuracy: 0.271844660194
k=100
Accuracy: 0.26213592233
k=150
Accuracy: 0.252427184466
k=200
Accuracy: 0.242718446602
k=400
Accuracy: 0.242718446602

50から100のあたりでピークがあります。

8.5.4 Modify the Bayes classifier to use a different probability model than Gaussian distributions. For example, try using the frequency counts of each feature in the training set. Compare the results to using a Gaussian distribution for some different datasets.

ガウス分以外の確率分布を試せ

私の回答

ポアッソン分布を試してみます
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%204-1.ipynb

from numpy import *
from PIL import *
import pickle
from pylab import *
import os
import pickle
import bayes
import imtools

with open('points_normal.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

bc = bayes.BayesClassifier()
bc.train([class_1, class_2], [1, -1])

with open('points_normal_test.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

print bc.classify(class_1[:10])[0]
print bc.classify(class_2[:10])[0]

# draw points and boundary line
def classify(x, y, bc=bc):
    points = vstack((x, y))
    return bc.classify(points.T)[0]

imtools.plot_2D_boundary([-6, 6, -6, 6], [class_1, class_2], classify, [1, -1])
show()

Gaussian
[1 1 1 1 1 1 1 1 1 1]
[-1 -1 -1 -1 -1 -1 -1 -1 -1 -1]

with open('points_normal.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)
import heurisitc
bcp = poisson.BayesClassifier()
bcp.train([class_1, class_2], [1, -1])
with open('points_normal_test.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)
print bcp.classify(class_1[:10])[0]
print bcp.classify(class_2[:10])[0]
def classify(x, y, bc=bcp):
    points = vstack((x, y))
    return bc.classify(points.T)[0]

imtools.plot_2D_boundary([-6, 6, -6, 6], [class_1, class_2], classify, [1, -1])
show()

Poisson
[1 1 1 1 1 1 1 1 1 1]
[-1 -1 -1 -1 -1 -1 -1 -1 -1 -1]

with open('points_ring.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

bc = bayes.BayesClassifier()
bc.train([class_1, class_2], [1, -1])

with open('points_ring_test.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

def classify(x, y, bc=bc):
    points = vstack((x, y))
    return bc.classify(points.T)[0]

imtools.plot_2D_boundary([-6, 6, -6, 6], [class_1, class_2], classify, [1, -1])
show()

Gaussian

import poisson
poisson = reload(poisson)

with open('points_ring.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

bcp = poisson.BayesClassifier()
bcp.train([class_1, class_2], [1, -1])

with open('points_ring_test.pkl', 'r') as f:
    class_1 = pickle.load(f)
    class_2 = pickle.load(f)
    labels = pickle.load(f)

def classify(x, y, bc=bcp):
    points = vstack((x, y))
    return bc.classify(points.T)[0]

imtools.plot_2D_boundary([-6, 6, -6, 6], [class_1, class_2], classify, [1, -1])
show()

Poisson

ポアッソン分布でできた境界が直線なのが意外です。あってるんでしょうか。

8.5.5 Experiment with using non-linear SVMs for the gesture recognition problem. Try polynomial kernels and increase the degree (using the "-d" parameter) incrementally. What happens to the classification performance on the training set and the test set. With a non-linear classifier there is a risk of training and optimizing it for a specific set so that performance is close to perfect on the training set but the classifier has poor performance on other test sets. This phenomenon of breaking the generalization capabilities of a classifier is called overfitting and should be avoided.

手話認識を非線形SVMでやりなおせ

私の回答

ヒントとして多項式カーネルを指定して次数を順番に増やしてみよ、とあるのでそのようにします。
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%205.ipynb

これも内容は教科書の内容と似ているのでコードは割愛して結果だけ。
次数3から12まで
Accuracy: 0.793269230769
Accuracy: 0.793269230769
Accuracy: 0.798076923077
Accuracy: 0.798076923077
Accuracy: 0.793269230769
Accuracy: 0.783653846154
Accuracy: 0.778846153846
Accuracy: 0.769230769231
Accuracy: 0.764423076923

グラフ

５、６次にピークがあり、あとは落ちていきます。

8.5.6 Try some more advanced feature vectors for the sudoku character recognition problem. If you need inspiration, look at [4].

より高度な方法で数独の文字認識を行え

私の回答

https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%206-1.ipynb
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%206-2.ipynb

来ました、ざっくりした問題！
普通の問題ばかりで油断していました。いったい「高度な」というのはどの程度をさすんでしょうか？解釈次第によっては一大プロジェクトになりかねません。今回はただの演習問題としてやるので、ヒントのとおり参考文献からパクらせていただきます。

論文の内容にしたがって、ノイズ除去や閾値処理やサイズの調整を行った上で、黒字の割合、列方向の合計の最大値（各列の画素値の合計のうち最大値）、行方向の合計（各行の画素値の合計のうち最大値）、縦・横それぞれ白黒反転の回数、穴の数、モーメント、対称性、を計算し特徴値とします。

def compute_feature(im):
    """ Returns a feature vector for an
    ocr image patch. """

    im = imresize(im, (300, 300))
    im = im[30:-30, 30:-30]
    im = 1*(im<128)
    im_open = morphology.binary_opening(im, ones((9, 9)), iterations=2)
    im_open = 1*(im_open>0.5)

    col = np.sum(im_open, axis=0)
    row = np.sum(im_open, axis=1)
    try:
        xmin = max(0, col.nonzero()[0][0]-10)
        xmax = min(299, col.nonzero()[0][-1]+10)
        ymin = max(0, row.nonzero()[0][0]-10)
        ymax = min(299, row.nonzero()[0][-1]+10)
    except:
        xmin = 0
        xmax = 299
        ymin = 0
        ymax = 299

    im_crop = im_open[ymin:ymax, xmin:xymax]
    im_crop = imresize(im_crop, (24, 24))
    im_crop = 1*(im_crop>128)

    black_proportion = 1.0*sum(im_crop)/(im_crop.shape[0]*im_crop.shape[1])

    col_sum = np.sum(im_crop, axis=0)
    row_sum = np.sum(im_crop, axis=1)

    max_col_sum = 1.0*max(col_sum)/im_crop.shape[0]
    max_row_sum = 1.0*max(row_sum)/im_crop.shape[1]

    nbr_alt_cols = 1.0*count_alternative_lines(col_sum, len(row_sum)/2)/len(col_sum)
    nbr_alt_rows = 1.0*count_alternative_lines(row_sum, len(col_sum)/2)/len(row_sum)

    labels, nbr_of_objects = measurements.label(1-im_crop)
    nbr_of_holes = nbr_of_objects-1

    moments = get_moments(im_crop)

    symmetry = get_symmetry(im_crop)

    return [black_proportion, max_col_sum, max_row_sum, nbr_alt_cols, nbr_alt_rows,
           nbr_of_holes, moments[0], moments[1], moments[2], moments[3], moments[4],
           symmetry[0], symmetry[1]]

この特徴値で文字判定すると、
Accuracy = 77.7332% (775/997) (classification)
と、なりました。あまり高くはありませんが、１０の特徴値しか使っていないと考えるとわるくありません。

実際に数独に適用してみましょう

imname = 'sudoku_images/sudokus/sudoku18.JPG'
vername = 'sudoku_images/sudokus/sudoku18.sud'
im = array(Image.open(imname).convert('L'))

x = find_sudoku_edges(im, axis=0)
y = find_sudoku_edges(im, axis=1)

crops = []
for col in range(9):
    for row in range(9):
        crop = im[y[col]:y[col+1], x[row]:x[row+1]]
        crops.append(compute_feature(crop))

res = svm_predict(loadtxt(vername), map(list, crops), m)[0]
res_im = array(res).reshape(9, 9)

figure()
gray()
imshow(im)
axis('off')
show()

print 'Result:'
print res_im

Accuracy = 100% (81/81) (classification)

Result:
[[ 0. 0. 0. 0. 1. 7. 0. 5. 0.]
[ 9. 0. 3. 0. 0. 5. 2. 0. 7.]
[ 0. 0. 0. 0. 0. 0. 4. 0. 0.]
[ 0. 1. 6. 0. 0. 4. 0. 0. 2.]
[ 0. 0. 0. 8. 0. 1. 0. 0. 0.]
[ 8. 0. 0. 5. 0. 0. 6. 4. 0.]
[ 0. 0. 9. 0. 0. 0. 0. 0. 0.]
[ 7. 0. 2. 1. 0. 0. 8. 0. 9.]
[ 0. 5. 0. 2. 3. 0. 0. 0. 0.]]

いい感じです。

8.5.7 . Implement a method for automatically aligning the sudoku grid. Try for example feature detection with RANSAC, line detection or detecting the cells using morphological and measurement operations from scipy.ndimage (http://docs.scipy.org/doc/scipy/reference/ndimage.html). Bonus task, solve the rotation ambiguity of finding the "up" direction. For example, you could try rotating the rectified grid and let the OCR classifier’s accuracy vote for the best orientation.

枠線を自動的に位置合わせせよ

私の回答

これもざっくりしていますが、やってみましょう。
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%207.ipynb

ソースコードは長過ぎるのでgithubを参照下さい。以下、やったことの解説と結果。

この図(sudoku7.JPG)を例として使います

まずエッジ強調して閾値をつかって２値化します

次にハフ変換して、直線のパラメーターを推定します。

直線を当てはめてみましょう

交点を見つけます

四隅を見つけます

この情報を元にndimageのgeometric_transformに渡すパラメータを作り、変形します。

上手くいきました。

次に回転の曖昧さを取り除くために、９０度ずつ回転させてOCRにかけ、もっとも精度の高そうなところを上方向に採用します。これは問題分のヒントのとおりです。

Accuracy = 92.5926% (75/81) (classification)
Accuracy = 40.7407% (33/81) (classification)
Accuracy = 71.6049% (58/81) (classification)
Accuracy = 43.2099% (34/81) (classification)
Best accuracy: 0.925925925926
[[ 0. 0. 0. 0. 0. 5. 0. 3. 0.]
[ 3. 0. 3. 0. 0. 4. 0. 4. 0.]
[ 0. 3. 0. 0. 3. 0. 2. 1. 0.]
[ 0. 0. 7. 0. 6. 2. 0. 9. 0.]
[ 0. 0. 0. 8. 0. 9. 0. 0. 0.]
[ 0. 6. 0. 7. 5. 0. 3. 0. 0.]
[ 0. 9. 6. 0. 4. 0. 0. 3. 0.]
[ 0. 2. 0. 6. 0. 0. 5. 0. 8.]
[ 0. 3. 0. 2. 0. 0. 0. 0. 0.]]

うまく行きました。
あれ？でも精度を知るにはどんな字があるかわからないといけませんね。このデータベースは文字があらかじめわかっていますが、それは実際の応用では非現実的な仮定です。
うーん、まだちょっと足りなかったようです。

8.5.8 For a more challenging classification problem than the sudoku digits, take a look at the MNIST database of handwritten digits http://yann.lecun.com/exdb/ mnist/. Try to extract some features and apply SVM to that set. Check where your performance ends up on the ranking of best methods (some are insanely good).

MINSTの手書き文字認識を試してみよ

私の回答

最近ではMINSTも時代遅れなようですが、演習問題なのでそのあたりは問題ありません。
https://github.com/moizumi99/CVBookExercise/blob/master/Chapter-8/CV%20Book%20Chapter%208%20Exercise%208.ipynb

やってみました。

from numpy import *
from PIL import *
import pickle
from pylab import *
import os
from scipy.misc import *
from matplotlib.pyplot import *
# data is taken from this page
# http://yann.lecun.com/exdb/mnist/
train_file = list(fromfile('handwriting/train-images.idx3-ubyte', uint8, -1))

nbr_img = train_file[4]*(256**3)+train_file[5]*(256**2)+train_file[6]*256+train_file[7]
width = train_file[8]*(256**3)+train_file[9]*(256**2)+train_file[10]*256+train_file[11]
height = train_file[12]*(256**3)+train_file[13]*(256**2)+train_file[14]*256+train_file[15]

features = []
siz = height*width
for i in range(nbr_img):
    img = train_file[16+i*siz:16+(i+1)*siz]
    features.append(img)

label_file = list(fromfile('handwriting/train-labels.idx1-ubyte', uint8, -1))
nbr_labels = label_file[4]*(256**3)+label_file[5]*(256**2)+label_file[6]*256+label_file[7]
labels = label_file[8:8+nbr_labels]

test_file = list(fromfile('handwriting/t10k-images.idx3-ubyte', uint8, -1))
nbr_test_img = test_file[4]*(256**3)+test_file[5]*(256**2)+test_file[6]*256+test_file[7]
width = test_file[8]*(256**3)+test_file[9]*(256**2)+test_file[10]*256+test_file[11]
height = test_file[12]*(256**3)+test_file[13]*(256**2)+test_file[14]*256+test_file[15]
test_features = []
siz = height*width
for i in range(nbr_test_img):
    img = test_file[16+i*siz:16+(i+1)*siz]
    test_features.append(img)

test_label_file = list(fromfile('handwriting/t10k-labels.idx1-ubyte', uint8, -1))
nbr_test_labels = test_label_file[4]*(256**3)+test_label_file[5]*(256**2)+test_label_file[6]*256+test_label_file[7]
test_labels = test_label_file[8:8+nbr_test_labels]

from svmutil import *

nbr_test = 30000

prob = svm_problem(labels[:nbr_test], features[:nbr_test])
param = svm_parameter('-t 0')
m = svm_train(prob, param)
res = svm_predict(labels[:nbr_test], features[:nbr_test], m)
nbr_test2 = 10000
res = svm_predict(test_labels[:nbr_test2], test_features[:nbr_test2], m)

Accuracy = 91.43% (9143/10000) (classification)

結果91%ほどですが、ランキングの上位よりはだいぶ低いようです。

8.5.9 If you want to dive deeper in classifiers and machine learning algorithms, take a look at the scikit.learn package (http://scikit-learn.org/) and try some of the algorithms on the data in this chapter.

scikit.learnパッケージを調べ、アルゴリズムをいくつか試してみろ

私の回答

試してみましたが、あまり面白い結果はでませんでした。
したがって、結果は割愛します。

８章終わり

そろそろ終盤です。だんだんざっくりした問題も増えてきましたが、めげずにいきましょう。