人工知能に関する断創録

このブログでは人工知能のさまざまな分野について調査したことをまとめています(更新停止: 2019年12月31日)

Visual Wordsを用いた類似画像検索

に続くOpenCVプロジェクト第三弾です。今回は、上の二つをふまえてカラーヒストグラムではなく、局所特徴量(SIFTやSURF)を用いた類似画像検索を試してみます。局所特徴量はグレースケール画像から抽出するため、カラーヒストグラムと違って色は見ていません。画像の模様(テクスチャ)で類似性を判定します。

実験環境は、Windows 7、MinGW C++コンパイラ、OpenCV2.0、Python 2.5です。EclipseでMinGWを使う方法はEclipseでOpenCV(2009/10/16)を参照してください。Visual C++にはないディレクトリスキャン関数を一部使っているのでVisual C++を使う場合は、少しだけ修正が必要です。

Bag-of-Visual Words

まず、画像のBag-of-Visual Words表現について簡単にまとめておきます。Bag-of-Visual Wordsを用いると画像をただ1つの特徴ベクトルで表現できます。Bag-of-Visual Wordsは、別名としてBag-of-FeaturesBag-of-Keypointsという呼び方もあります。Googleで調べてみたところ

  • Bag-of-Features(1,660,000件)
  • Bag-of-Visual Words(113,000件)
  • Bag-of-Keypoints(58,400件)

となったのでBag-of-Featuresという呼び方が一番多いんでしょうか?私自身はBag-of-Visual Wordsという呼び方が一番好きです。この手法はもともと自然言語処理でよく使われるBag-of-Wordsのアナロジーとして画像へ応用した手法だからです。

Bag-of-Wordsは文書をベクトルで表現する手法で自然言語処理の分野では当たり前のように使います。簡単に言うと、文書を単語の集合として扱います。実際は、対象文書集合からTF-IDFで特徴語を抽出して文書ベクトルの次元とし、TF-IDFの重みを文書ベクトルの値とする方法がよく使われます。あとでナイーブベイズを利用した文書分類を紹介する予定ですがそのときに使ってみます。

f:id:aidiary:20100227224223p:plain

では、画像における単語とは何か?という話になるのですが、それは今まで使ってきた局所特徴量です。Bag-of-Visual Wordsは、テキストの単語⇔画像の局所特徴量と対応させる方法です。上図のようにまずVisual Wordsを構築する画像集合を用意します。そして、それぞれの画像から局所特徴量を抽出します。図では、わかりやすいように局所特徴量を小さな画像で表現していますが、実際は128次元のSIFTやSURFの特徴ベクトルです。そして、全部の特徴ベクトルをK個のクラスタにクラスタリングします。クラスタリングにはK-meansという手法がよく使われます。そして、K個のクラスタの各セントロイド(中心となるベクトル)をそれぞれVisual Wordとします。セントロイドも特徴ベクトルと同じく128次元ベクトルで表されます。つまり、K個のクラスタにクラスタリングしたらK個のVisual Wordsができます。この例では、K=7なので7個のVisual Wordが得られます。

f:id:aidiary:20100227224224p:plain

次に、上図のように画像データベース中の全画像(Visual Wordsを求めたものも含む)をVisual Wordsを次元とするヒストグラムに変換します。これは、画像中の各局所特徴量について一番近いVisual Wordを検索し、そのVisual Wordに投票することで行います。最終的にこのヒストグラムが画像の特徴ベクトルになります。この例では、7次元の特徴ベクトルです。実際は、Visual Wordは数百から数千個選ぶので数百から数千次元のベクトルになります。Visual Wordsを使うと一枚の画像を数百次元のベクトル一本で表せます。画像を一本のベクトルで表せると、Support Vector Machineなどの分類器が使いやすくなってうれしいですね。このヒストグラムの形は類似画像ほど似ているはずです。たとえば、上図だとアコーディオン同士のヒストグラムは似ていますね。ヒストグラム同士の類似度は、類似画像検索システムを作ろう(2009/10/3)で使ったヒストグラムインターセクションなどが使えます。

説明が長くなりましたが、実際に試してみます。もう少し詳しい説明は参考文献を参照してください。

画像データセットを用意

実験で使う画像データセットを用意します。今回は類似画像検索システムを作ろう(2009/10/3)で使ったCaltech101をまた使います。ダウンロードしたデータを解凍すると101_ObjectCategoriesというフォルダの中に101個のカテゴリフォルダ(102個あるけどBACKGROUND_Googleが余計?)があり、その中に画像がたくさんあります。たとえば、airplanes フォルダ(カテゴリ)には飛行機の写真がたくさん入ってます。

まず、BACKGROUND_Google, Faces, Faces_easyは少し特殊な画像なのでフォルダごと削除します。次に、フォルダごとに画像ファイルが分散されてるので一括処理しやすいようにcaltech101というフォルダに画像ファイルだけまとめます。また各カテゴリフォルダ内の画像ファイル名はimage_0001.jpgとかついてます。airplanes(飛行機)にもdolphin(イルカ)にも同じファイル名がついてるので区別できなくて面倒です。そこで、下のようにカテゴリ名-数字.jpgというファイル名に変換します。

    airplanes/image_0001.jpg -> caltech101/airplanes-0001.jpg
    dolphin/image_0001.jpg -> caltech101/dolphin-0001.jpg

一括して変換するpythonスクリプトです。

#coding:utf-8
import codecs
import os
import os.path
import shutil

TARGET = "101_ObjectCategories"
OUTDIR = "caltech101"

# 出力ディレクトリがなければ作る
if not os.path.exists(OUTDIR):
    os.mkdir(OUTDIR)

for category in os.listdir(TARGET):
    for file in os.listdir("%s/%s" % (TARGET, category)):
        # 101_ObjectCategories/airplanes/image_0001.jpg
        image_file = "%s/%s/%s" % (TARGET, category, file)
        # caltech101/airplanes-0001.jpg
        rename_file = "%s/%s-%s" % (OUTDIR, category, file.replace("image_", ""))
        # ファイルをコピー
        print "%s -> %s" % (image_file, rename_file)
        shutil.copyfile(image_file, rename_file)

これで、9145枚の画像ファイルだけcaltech101フォルダにリネイムコピーされます。これでヒストグラムの抽出処理などやりやすくなります。次に今回の実験は非常に計算量が大きいのでaccordion, bonsai, cougar_face, dalmatian, dollar_bill, euphonium, hedgehog, grand_piano, Motorbikes, yin_yangの10カテゴリの1番目から50番目の画像のみ抽出してcaltech10という500枚の画像からなるサブセットデータ集合を作ります。

#coding:utf-8
import os.path
import shutil

IMAGE_DIR = "caltech101"
SUBSET_DIR = "caltech10"

# 対象カテゴリ
TARGET_CATEGORY = ["accordion", "bonsai", "cougar_face", "dalmatian", "dollar_bill",
                   "euphonium", "hedgehog", "grand_piano", "Motorbikes", "yin_yang", ]

# ファイル番号が50以下の画像のみ対象
TOP = 50

# 出力ディレクトリがなければ作る
if not os.path.exists(SUBSET_DIR):
    os.mkdir(SUBSET_DIR)

for file in os.listdir(IMAGE_DIR):
    try:
        cat = file.split("-")[0]            # ファイルのカテゴリ名を取得
        num = int(file.split("-")[1][0:4])  # ファイル番号を取得
    except:
        continue
    
    # 対象カテゴリで数字がTOP以下のファイルのみコピー
    if cat in TARGET_CATEGORY and num <= TOP:  # 対象カテゴリの場合
        source_image = "%s/%s" % (IMAGE_DIR, file)
        dest_image = "%s/%s" % (SUBSET_DIR, file)
        shutil.copyfile(source_image, dest_image)

全体の流れ

画像データセットが用意できたんで、C++とOpenCVを使ってVisual Wordsを用いた類似画像検索を実装していきます。全体的な流れをまとめると

  1. 全画像からSURF局所特徴量を抽出
  2. 全画像のSURFベクトルをK-meansでK個のクラスタにクラスタリングして各クラスタのセントロイド(Visual Word)を求める
  3. 全画像をヒストグラムに変換する

です。ヒストグラムは後で使えるようにファイルに出力します。

visual_words.cpp

まず、メイン関数から。

#include <cv.h>
#include <highgui.h>
#include <iostream>
#include <fstream>
#include <dirent.h>

using namespace std;

const char* IMAGE_DIR = "caltech10";
const int DIM = 128;
const int SURF_PARAM = 400;
const int MAX_CLUSTER = 500;  // クラスタ数 = Visual Wordsの次元数

int main() {
    int ret;

    // IMAGE_DIRの各画像から局所特徴量を抽出
    cout << "Load Descriptors ..." << endl;
    CvMat samples;
    vector<float> data;
    ret = loadDescriptors(samples, data);

    // 局所特徴量をクラスタリングして各クラスタのセントロイドを計算
    cout << "Clustering ..." << endl;
    CvMat* labels = cvCreateMat(samples.rows, 1, CV_32S);        // 各サンプル点が割り当てられたクラスタのラベル
    CvMat* centroids = cvCreateMat(MAX_CLUSTER, DIM, CV_32FC1);  // 各クラスタの中心(セントロイド) DIM次元ベクトル
    cvKMeans2(&samples, MAX_CLUSTER, labels, cvTermCriteria(CV_TERMCRIT_EPS+CV_TERMCRIT_ITER, 10, 1.0), 1, 0, 0, centroids, 0);
    cvReleaseMat(&labels);  // ラベルは使わない

    // 各画像をVisual Wordsのヒストグラムに変換する
    // 各クラスターの中心ベクトル、centroidsがそれぞれVisual Wordsになる
    cout << "Calc Histograms ..." << endl;
    calcHistograms(centroids);
    cvReleaseMat(&centroids);

    return 0;
}

(1)局所特徴量の抽出と(3)ヒストグラムの抽出は長いのでそれぞれloadDescriptors()、calcHistograms()という関数にまとめてあります。クラスタリングは関数を作るほど長くないのでそのままです。K-meansクラスタリングは、cvKMeans2()というOpenCVの関数で簡単にできます。samplesは全画像の局所特徴量の行列(特徴ベクトル数x128次元)です。MAX_CLUSTERはクラスタ数で適当に500としています。つまり、Visual Wordsは500個できます。labelsは今回は使いませんが、各特徴ベクトルがどのクラスタに割り当てられたかを表します。centroidは各クラスタのセントロイドベクトルが格納される500クラスタx128次元の行列です。

局所特徴量を抽出

caltech10の500枚の画像から局所特徴量SURFを抽出します。SURFの抽出方法は、SURFの抽出(2009/10/30)を参照してください。

/**
 * 画像ファイルからSURF特徴量を抽出する
 * @param[in]  filename            画像ファイル名
 * @param[out] imageKeypoints      キーポイント
 * @param[out] imageDescriptors    各キーポイントのSURF特徴量
 * @param[out] storage             メモリ領域
 * @return 成功なら0、失敗なら1
 */
int extractSURF(const char* filename, CvSeq** keypoints, CvSeq** descriptors, CvMemStorage** storage) {
    // グレースケールで画像をロードする
    IplImage* img = cvLoadImage(filename, CV_LOAD_IMAGE_GRAYSCALE);
    if (img == NULL) {
        cerr << "cannot load image: " << filename << endl;
        return 1;
    }

    *storage = cvCreateMemStorage(0);
    CvSURFParams params = cvSURFParams(SURF_PARAM, 1);
    cvExtractSURF(img, 0, keypoints, descriptors, *storage, params);

    return 0;
}

/**
 * IMAGE_DIRにある全画像から局所特徴量を抽出し行列へ格納する
 * @param[out]   samples    局所特徴量の行列
 * @param[out]   data       samplesのデータ領域
 * @return 成功なら0、失敗なら1
 */
int loadDescriptors(CvMat& samples, vector<float>& data) {
    // IMAGE_DIRを開く
    DIR* dp;
    if ((dp = opendir(IMAGE_DIR)) == NULL) {
        cerr << "cannot open directory: " << IMAGE_DIR << endl;
        return 1;
    }

    // IMAGE_DIRの画像ファイル名を走査
    struct dirent* entry;
    while ((entry = readdir(dp)) != NULL) {
        char* filename = entry->d_name;
        if (strcmp(filename, ".") == 0 || strcmp(filename, "..") == 0) {
            continue;
        }

        // パス名に変換
        // XXX.jpg -> IMAGE_DIR/XXX.jpg
        char filepath[256];
        snprintf(filepath, sizeof filepath, "%s/%s", IMAGE_DIR, filename);

        // SURFを抽出
        CvSeq* keypoints = NULL;
        CvSeq* descriptors = NULL;
        CvMemStorage* storage = NULL;
        int ret = extractSURF(filepath, &keypoints, &descriptors, &storage);
        if (ret != 0) {
            cerr << "error in extractSURF" << endl;
            return 1;
        }

        // ファイル名と局所特徴点の数を表示
        cout << filepath << "\t" << descriptors->total << endl;

        // 特徴量を構造化せずにdataへ追加
        for (int i = 0; i < descriptors->total; i++) {
            float* d = (float*)cvGetSeqElem(descriptors, i);  // 128次元ベクトル
            for (int j = 0; j < DIM; j++) {
                data.push_back(d[j]);
            }
        }

        cvClearSeq(keypoints);
        cvClearSeq(descriptors);
        cvReleaseMemStorage(&storage);
    }

    // dataをCvMat形式に変換
    // CvMatはdataを参照するためdataは解放されないので注意
    int rows = data.size() / DIM;  // CvMatの行数(=DIM次元特徴ベクトルの本数)
    cvInitMatHeader(&samples, rows, DIM, CV_32FC1, &data[0]);

    return 0;
}

loadDescriptors()でディレクトリスキャンの関数を使っていますが、これはVisual C++にはないと思います。似たような関数があるはずなので修正してください。

画像をヒストグラムに変換

次に、500枚の画像をVisual Wordsを次元としたヒストグラムに変換する関数です。ヒストグラムはタブ区切り形式でhistograms.txtに出力されます。画像の各局所特徴量にもっとも類似した(最近傍の)Visual Wordを見つけると書きましたが、Visual Wordsをkd-treeでインデキシングして検索速度を上げています。kd-treeを用いた高速化は、最近傍探索の高速化(2009/12/12)を参照してください。

/**
 * IMAEG_DIRの全画像をヒストグラムに変換して出力
 * 各画像の各局所特徴量を一番近いVisual Wordsに投票してヒストグラムを作成する
 * @param[in]   visualWords     Visual Words
 * @return 成功なら0、失敗なら1
 */
int calcHistograms(CvMat* visualWords) {
    // 一番近いVisual Wordsを高速検索できるようにvisualWordsをkd-treeでインデキシング
    CvFeatureTree* ft = cvCreateKDTree(visualWords);

    // 各画像のヒストグラムを出力するファイルを開く
    fstream fout;
    fout.open("histograms.txt", ios::out);
    if (!fout.is_open()) {
        cerr << "cannot open file: histograms.txt" << endl;
        return 1;
    }

    // IMAGE_DIRの各画像をヒストグラムに変換
    DIR* dp;
    if ((dp = opendir(IMAGE_DIR)) == NULL) {
        cerr << "cannot open directory: " << IMAGE_DIR << endl;
        return 1;
    }

    struct dirent* entry;
    while ((entry = readdir(dp)) != NULL) {
        char* filename = entry->d_name;
        if (strcmp(filename, ".") == 0 || strcmp(filename, "..") == 0) {
            continue;
        }

        char filepath[256];
        snprintf(filepath, sizeof filepath, "%s/%s", IMAGE_DIR, filename);

        // ヒストグラムを初期化
        int* histogram = new int[visualWords->rows];
        for (int i = 0; i < visualWords->rows; i++) {
            histogram[i] = 0;
        }

        // SURFを抽出
        CvSeq* keypoints = NULL;
        CvSeq* descriptors = NULL;
        CvMemStorage* storage = NULL;
        int ret = extractSURF(filepath, &keypoints, &descriptors, &storage);
        if (ret != 0) {
            cerr << "error in extractSURF" << endl;
            return 1;
        }

        // kd-treeで高速検索できるように特徴ベクトルをCvMatに展開
        CvMat* mat = cvCreateMat(descriptors->total, DIM, CV_32FC1);
        CvSeqReader reader;
        float* ptr = mat->data.fl;
        cvStartReadSeq(descriptors, &reader);
        for (int i = 0; i < descriptors->total; i++) {
            float* desc = (float*)reader.ptr;
            CV_NEXT_SEQ_ELEM(reader.seq->elem_size, reader);
            memcpy(ptr, desc, DIM*sizeof(float));
            ptr += DIM;
        }

        // 各局所特徴点についてもっとも類似したVisual Wordsを見つけて投票
        int k = 1;  // 1-NN
        CvMat* indices = cvCreateMat(keypoints->total, k, CV_32SC1);  // もっとも近いVisual Wordsのインデックス
        CvMat* dists = cvCreateMat(keypoints->total, k, CV_64FC1);    // その距離
        cvFindFeatures(ft, mat, indices, dists, k, 250);
        for (int i = 0; i < indices->rows; i++) {
            int idx = CV_MAT_ELEM(*indices, int, i, 0);
            histogram[idx] += 1;
        }

        // ヒストグラムをファイルに出力
        fout << filepath << "\t";
        for (int i = 0; i < visualWords->rows; i++) {
            fout << float(histogram[i]) / float(descriptors->total) << "\t";
        }
        fout << endl;

        // 後始末
        delete[] histogram;
        cvClearSeq(keypoints);
        cvClearSeq(descriptors);
        cvReleaseMemStorage(&storage);
        cvReleaseMat(&mat);
        cvReleaseMat(&indices);
        cvReleaseMat(&dists);
    }

    fout.close();
    cvReleaseFeatureTree(ft);

    return 0;
}

作成されたヒストグラムファイルをはっときます。
histograms.txt

Visual Wordsのヒストグラムを用いた類似画像検索

では、histograms.txtを使って類似画像検索のプログラムを書いてみます。もうOpenCVは使わないのでくそ面倒なC++を使う必要ないですね(笑)Pythonでいきます。

画像の類似度にはHistogram Intersectionを使います。Histogram Intersectionは、類似画像検索システムを作ろう(2009/10/3)を参照してください。次のスクリプトは、クエリとして与えた画像とその他全部の画像とのHistogram Intersectionを計算し、類似度が高い順に上位10位を返すPythonスクリプトです。Pythonの画像処理ライブラリPIL (Python Image Library)を使うので別途インストールしてください。全画像のヒストグラムはあらかじめメモリにロードしておきます

#coding:utf-8
import codecs
import os
import sys
from PIL import Image

IMAGE_DIR = "caltech10"
HIST_DIR = "hist"

# ヒストグラムをロードする
def load_hist():
    hist = {}
    fp = open("histograms.txt", "r")
    for line in fp:
        line = line.rstrip()
        data = line.split("\t")
        file = data[0]
        h = [float(x) for x in data[1:]]
        hist[file] = h
    fp.close()
    return hist

# 正規化Histogram Intersectionを計算する
def calc_hist_intersection(hist1, hist2):
    total = 0
    for i in range(len(hist1)):
        total += min(hist1[i], hist2[i])
    return float(total) / sum(hist1)

def main():
    hist = load_hist()

    while True:
        # クエリとなるヒストグラムファイル名を入力
        query_file = raw_input("query? > ")
        
        # 終了
        if query_file == "quit":
            break
        
        # パスに変換
        query_file = IMAGE_DIR + "/" + query_file
        
        # 存在しないヒストグラムファイル名のときは戻る
        if not hist.has_key(query_file):
            print "no histogram"
            continue
        
        # クエリと他の全画像の間で類似度を計算
        result = []
        query_hist = hist[query_file]
        for target_file in hist.keys():
            target_hist = hist[target_file]
            d = calc_hist_intersection(query_hist, target_hist)
            result.append((d, target_file))
        
        # 類似度の大きい順にソート
        result.sort(reverse=True)
        
        # 上位10位を表示(1位はクエリ画像)
        # PILを使って300x300を1単位として2行5列で10個の画像を並べて描画
        p = 0
        canvas = Image.new("RGB", (1500, 600), (255,255,255))  # 白いキャンバス
        for score, filename in result[0:10]:
            print "%f\t%s" % (score, filename)
            img = Image.open(filename)
            pos = (300*(p%5), 300*(p/5))
            canvas.paste(img, pos)
            p += 1
        canvas.resize((1500/2, 600/2))
        canvas.show()
        canvas.save("result.jpg", "JPEG")

if __name__ == "__main__":
    main()

実行すると、query?と聞いてくるのでcaltech10フォルダの画像ファイル名をクエリとして与えます。この結果としてクエリと類似する上位10個の画像が表示され、result.jpgにも保存します。

query? > accordion-0001.jpg
1.000000	caltech10/accordion-0001.jpg
0.762976	caltech10/accordion-0049.jpg
0.757379	caltech10/accordion-0025.jpg
0.753379	caltech10/accordion-0036.jpg
0.748672	caltech10/accordion-0042.jpg
0.739771	caltech10/accordion-0022.jpg
0.737838	caltech10/accordion-0002.jpg
0.734934	caltech10/cougar_face-0044.jpg
0.733918	caltech10/accordion-0037.jpg
0.733821	caltech10/accordion-0010.jpg

実験結果

では、いくつかのクエリ画像を与えて試してみます。左上にあるのがクエリ画像です。

  • アコーディオン

f:id:aidiary:20100227210336j:plain

おー、1個変なの入っているけどちゃんと検索できてますね。特定物体認識と違って、完全に同じ画像でなくても(特定物体認識でも回転、スケール変化などには対応できますが)ちゃんと類似検索できるってのがVisual Wordsを使う利点ですね。じゃ、別の例もいくつか。

  • 盆栽

f:id:aidiary:20100227210337j:plain

  • クーガー

f:id:aidiary:20100227210338j:plain

  • ダルメシアン

f:id:aidiary:20100227210339j:plain

  • ドル札

f:id:aidiary:20100227210340j:plain

  • ユーフォニアム

f:id:aidiary:20100227210341j:plain

  • グランドピアノ

f:id:aidiary:20100227210342j:plain

  • はりねずみ

f:id:aidiary:20100227210343j:plain

  • バイク

f:id:aidiary:20100227210344j:plain

  • 陰陽

f:id:aidiary:20100227210345j:plain

陰陽の結果を見ればわかるようにVisual Wordsは色の類似性は見ていないので注意してください。クーガーとかダルメシアンで漫画っぽい絵でも類似判定できてるのはすごいです。ドル札も額面は違うのにちゃんと類似判定できてます。はりねずみの中の陰陽もなんとなくわかりますね(笑)

(注)Visual Wordsの学習に使った画像とテスト用の画像が同じのは評価としてはちょっとまずいと思います。ちゃんと評価したい人は学習用画像セットとテスト用画像セットは分けましょう。

まとめ

今回は、局所特徴量を用いた類似画像検索を作りました。

上の結果を見ると何かすごいぞと思うのですが、実際はまったくダメダメなクエリ画像もありました。どういうときにうまくいかないかはもう少し調査が必要です。あとVisual Wordsの学習で使わなかったカテゴリ、たとえばイルカのクエリ画像を与えても類似画像は探せないと思います。Visual Wordsにイルカの特徴が含まれておらずヒストグラムがうまく作れないはずです(まだ試してないのでやってみてください)。

そんなわけでVisual Wordsで汎用的な類似画像検索システムを構成するのは難しいのではないか?と推測しています。実際に、Visual Wordsを類似画像検索に応用する研究は少ないように思います。Visual Wordsは主に一般物体認識で使われているようです。クエリとして与えた画像のカテゴリ(アコーディオンとかはりねずみとか)を分類するタスクです。Support Vector Machineを適用した研究をよく見かけます。いずれ一般物体認識も試してみたいですね。

Visual Wordsは、局所特徴量のクラスタリングが必要なわけですが、はっきり言ってK-meansがネックです。K-meansはすごい遅いです。bayonを使って画像からbag-of-keypointsを求める(この人のブログは愛読してます)を参考にCLUTObayonも検討していたのですが、いまいちデータの作り方がわかりませんでした・・・bayonは文書単語行列のようなスパース行列用のクラスタリングアルゴリズムだと理解したのですが、局所特徴量行列のような密な実数値ベクトルにも適用できるのでしょうか?bayonはもう少し詳しく調査してみる予定です。

参考文献