量子 SVM (QSVM) — 見慣れたデータセットを題材に

目的

この記事では量子 SVM について見慣れたものを用いた簡単な説明を行いたいと思います。

少し前に Qiskit の QSVM のコンテンツ Quantum feature maps and kernels（以降、文献 [QISKIT]）に取り組んでいて「難しいな」と思って苦しみました。そういう私のような普通のチュートリアルだと難しいと感じる人に「QSVM はそんなに怖くないよ」程度のことを伝えるのがゴールとなります。敬愛する加藤敏夫先生の著書「位相解析」のお言葉を拝借しますと「本記事のようなだらだらしたものが一つくらいあってもよいのではないかと考える」のです。

やること

• 文献 [QISKIT] を更に簡単に砕いて、同記事を読むときに「うっ・・・」とならない準備をする。（英語が「うっ・・・」となる場合、同コンテンツは邦訳済みなので言語切り替えをしてください・・・）
• 「量子コンピュータさっぱりわからん」人（より具体的には私自身）の心理的安全性を高める。
• 折角なので Blueqat SDK + scikit-learn だけで実装する。

やらないこと

• 量子アルゴリズムちょっとわかる人向けにスマートな説明をするということはしない。
• 発展的な最新のやり方は提供しない。
• カーネル法そのものの詳細な説明は行わない。

となります。また、基本的な道具として文献 [QISKIT] を大いに参考にします。

先に結論: 手短に「量子カーネル SVM」とは？

古典的なカーネル SVM において、グラム行列を用いた自作カーネルを適用したクラス分類問題の解き方があります。グラム行列はデータ間の内積をとるという古典計算で作成できますが、古典計算では近似しにくいような量子計算で求めた「グラム行列のようなもの」を自作カーネル（つまり「量子カーネル」）として用いて実行する SVM が量子カーネル SVM です。つまり、古典的なカーネル SVM との違いはカーネルの作り方が違うだけというのが私の認識です。———

量子カーネルを使うと何が嬉しいかは最後の「まとめ」の直前に書くことにします。

なおグラム行列は、データ \{\xi_i\}_{1 \leq i,j \leq N} がある時に、データ同士の内積を要素に持つ行列 G = (\xi_i \cdot \xi_j)_{1 \leq i,j \leq N} で定まる行列になります。ベクトルの内積というものは、仮にベクトルの長さを 1 だと思うことにするといわゆる「コサイン類似度」の式になります。よって、少々強引な言い方ですが「グラム行列はデータ間の類似度を成分に持つ行列である」となります。

SVM のおさらい

SVM とは、典型的には「犬」のデータと「猫」のデータのように 2 つのクラスからなるデータについてデータの特徴量を元にクラス分類を行う古典的な機械学習の手法の 1 つです。以下の図で赤い三角を「犬」、青い丸を「猫」とする時に、直線を引いて「犬」と「猫」の分類を行うことになります。

分類に先立って「犬」と「猫」をベクトル空間にマッピングしておく必要がありますが、そのためには何かしらの特徴、例えば

• 鼻の先から後頭部への長さ
• 胴体の長さに対する尻尾の長さの割合

といったものに基づいてプロットすることになります。

うまい特徴を選ぶことでこの図のように線を引いて区別できる場合は良いのですが、そうはいかない以下のようなケースもあります。

ちょっと難しいケース

これはよく見かけるデータセットですが、データセットの座標 (x,y) を (x,y,x^2+y^2) に変換してより高次元空間に埋め込み直すことで以下のように分類できるようになります。今度は空間なので平面で分類することになります。3 より大きい高次元空間での場合、超平面と呼んだりもします。

難しいケースでは具体的にどうするの？

元々データが存在していた空間を \Omega として、それをより高次元の空間 \mathcal{H} に写像してやると、データの分離が捗るのではないか？という考え方です。この同心円状のデータセットの場合、特徴写像 \Phi: (x,y) \mapsto (x,y,x^2+y^2) で 2 次元平面 \Omega から 3 次元空間 \mathcal{H} にうつしてあげました。

この手の設定では高次元の空間 \mathcal{H} は好きに選ばれるのではなく、“カーネル” と呼ばれるある写像 k(\cdot, \cdot) を持つようなものが選ばれます。\Omega のデータセットを \{\xi_i\}_{i=1}^N としますと、このカーネルというやつは

という関係を満たすことが知られています。ここで \langle \cdot, \cdot \rangle_\mathcal{H} は空間 \mathcal{H} の内積とします。重要なことは高次元の空間にデータを埋め込んだ時、埋め込んだデータ同士の内積はカーネルを通じて簡単に計算できるということです。

内積が簡単に計算できると何が嬉しいのかが見えないため、少しだけ補足します。それは、高次元空間でデータを分離する平面を求めるために以下のような (2) 式を扱う必要があり、式の最後のように高次元空間 \mathcal{H} での内積が絡んでくるのです。ちなみに式の細かい意味はどうでも良いです。

従いまして、内積の計算がコンピュータ上で高速に行えない場合、カーネル SVM は使い物にならないアルゴリズムということになるのですが、(1) 式のようにカーネルを通じて内積が高速に計算できるので嬉しいということになります。要約すると、カーネルを持つような空間を使うのは、技術的な都合です。この計算に都合のよい仕組みを提供する (1) 式の関係性のことはカーネルトリックなどと呼ばれます。

ちなみに、今回の場合、良いカーネルを持つ空間を見つけるのではなく、カーネル k を自分で作ってしまって「カーネル k に紐づく空間 \mathcal{H} と特徴写像 \Phi が存在する] と思うことにし、なおかつ「その \Phi は (1) 式を満たしているのだ」と考えることにします。つまり具体的な \mathcal{H} や \Phi は扱いません。もっとずるいことに、(2) 式よりデータセットを \Phi でうつした像の内積 \langle \Phi(\xi_i), \Phi(\xi_j) \rangle_\mathcal{H} さえ分かれば良いので、k という写像をつくる代わりに \boldsymbol{k(\xi_i, \xi_j)} の値だけを作成するということをします。

中間まとめ

• 元の空間 \Omega のデータセット \{\xi_i\}_{i=1}^N をクラス分類したい場合、元の空間のままでは難しい場合、高次元の空間 \mathcal{H} に埋め込めば解けるかもしれない。
• 高次元の空間 \mathcal{H} は “カーネル” と呼ばれる都合の良い写像に紐づけられたものを選ばないとコンピュータで低コストで計算できない。
• “カーネル” は自分で用意してしまって、具体的な \mathcal{H} は心の中だけにとどめてしまう。

だいぶ無理矢理ですが「カーネルを設計すれば良さそうだ」という雰囲気になってきたのではないでしょうか。

さて、以下では「鶏を割くに焉んぞ牛刀を用いん」になるのですが、この同心円のデータセットを量子カーネルを用いて分類します。「こんな簡単なデータセットに量子回路を使うのか！？」と思われるところですが、簡単なこともできないアルゴリズムに意味はないので、敢えて簡単なものを選びました。（初等的な解法としては、例えば文献 [PML] を参考に、多項式カーネルや RBF カーネルというありもののカーネルを用いれば簡単に解けます）

「難しいと感じる理論では、具体的で一番簡単なケースで理論の本質的な流れを追いかけろ」というのは私の恩師からの教えになります。

量子カーネルとは？

念押しになりましが、量子回路で計算したカーネルのことです。どういうものが良いのかはまだまだ研究中ということもあり、素直に論文 [H-C-T-H-K-C-G] のものを用います。Qiskit で描画したものを使うと、以下のような回路の断片を繋いだようなものになります。ここで H はアダマールゲート、P は位相ゲートになります。位相ゲートはグローバル位相を除いて R_z ゲートと同じですので、後の実装では R_z を用います。

ちなみにこの回路は Qisikit では ZZFeatureMap と呼ばれるものになるのですが、この名前はハミルトニアン Z \otimes Z の時間発展演算子 \exp (-i \theta Z \otimes Z) を含むからと考えられます。この時間発展のユニタリ演算子は R_{zz} ゲートとも呼ばれ、以下の形をしています。

この ZZFeatureMap にはパラメータの口 x[0] と x[1] が見えるのですがここにデータセット \{\xi_i\}_{i=0}^N の座標を埋め込みます。いわゆる「位相エンコーディング」と呼ばれるものです。具体的には \xi_i = (\xi_i^0, \xi_i^1) と 2 次元座標表示した場合に、x[0] のほうには \xi_i^0 を、x[1] のほうには \xi_i^1 を代入して使います。

そろそろ Python での実装にうつりたいと思います。まずはこの ZZFeatureMap を Blueqat SDK で実装します。なお、後で reps=2 として使うつもりですので、パラメータを部分適用した feature_map 関数もここで作っておきます。

import numpy as np
from functools import partial
from blueqat import Circuit

def zz_feature_map(x, reps):
    def sub_circuit(x):
        n_qubit = len(x)
        c = Circuit().h[:]
        for i in range(n_qubit):
            c.rz(2*x[i])[i]
        for i in range(n_qubit - 1):
            for j in range(i+1, n_qubit):
                c.cx[i, j].rz(2*(np.pi-x[i])*(np.pi-x[j]))[j].cx[i, j]
        return c

    c = Circuit()
    for _ in range(reps):
        c += sub_circuit(x)
    return c

feature_map = partial(zz_feature_map, reps=2)

折角なのでデータ \xi = (0, 0) を位相エンコーディングする形で回路を描画してみましょう。

c = zz_feature_map([0, 0], 1)
c.run(backend='draw')

<Figure size 1800x180 with 1 Axes>

どうでしょうか？Qiskit でいう ZZFeatureMap になっていると思います。先ほど注意しましたように、位相ゲート P と 回転ゲート R_z はグローバル位相を除いて同じゲートなので、ここでは R_z を用いています。

量子カーネルを作る

再び論文 [H-C-T-H-K-C-G]、特に p.5 を参考にして、上記の zz_feature_map を部品としてカーネルを作ります。この zz_feature_map を表す量子回路を \mathcal{U}_\Phi とします。論文によると、\Phi(\xi_i) = \mathcal{U}_\Phi (\xi_i) | 0^n \rangle と \Phi(\xi_j) = \mathcal{U}_\Phi (\xi_j) | 0^n \rangle の内積の絶対値の 2 乗

を計算して、これを (i,j)-成分に持つような行列

をカーネルとします。ここで n は量子ビット数ですが、今回はデータセットが 2 次元平面状にありましたので n=2 となります。 また、| 0^n \rangle = | 0 \rangle \otimes \cdots \otimes | 0 \rangle で、n 個の | 0 \rangle のテンソル積です。

冒頭で、量子計算で「グラム行列のようなもの」を自作カーネルとすると書きましたが、(4) 式がそれです。

簡単なキャッチフレーズを作りますと、「今回作ったカーネルはデータセットを高次元空間 \mathcal{H} に埋め込んだ時の類似度を成分に持つ行列である」となります。

さて、このカーネルを Blueqat SDK で実装しましょう。本来は \Phi(\xi_i) = \mathcal{U}_\Phi (\xi_i) の部分は測定を通じて確率振幅を求めるべきと思いますが、今回はずるをして状態ベクトルから計算します。つまり、numpy バックエンドを用いて .run() します。

def calculate_kernel(feature_map, x_data, y_data=None):
    if y_data is None:
        y_data = x_data
    x_matrix, y_matrix = [], []
    for x0, x1 in x_data:
        c = feature_map([x0, x1])
        sv = c.run(backend='numpy')
        x_matrix.append(sv)
    for y0, y1 in y_data:
        c = feature_map([y0, y1])
        sv = c.run(backend='numpy')
        y_matrix.append(sv)
    x_matrix, y_matrix = np.array(x_matrix), np.array(y_matrix)
    kernel = np.abs(
        y_matrix.conj() @ x_matrix.T
    )**2
    return kernel

今更ですが、データセットを用意します。

from sklearn.datasets import make_circles

def make_train_test_sets(test_ratio=.3):
    def split(arr, test_ratio):
        sep = int(arr.shape[0]*(1-test_ratio))
        return arr[:sep], arr[sep:]

    X, Y = make_circles(n_samples=200, noise=0.05, factor=0.4)

    A = X[np.where(Y==0)]
    B = X[np.where(Y==1)]

    A_label = np.zeros(A.shape[0], dtype=int)
    B_label = np.ones(B.shape[0], dtype=int)

    A_label = np.zeros(A.shape[0], dtype=int)
    B_label = np.ones(B.shape[0], dtype=int)
    A_train, A_test = split(A, test_ratio)
    B_train, B_test = split(B, test_ratio)
    A_train_label, A_test_label = split(A_label, test_ratio)
    B_train_label, B_test_label = split(B_label, test_ratio)
    X_train = np.concatenate([A_train, B_train])
    y_train = np.concatenate([A_train_label, B_train_label])
    X_test = np.concatenate([A_test, B_test])
    y_test = np.concatenate([A_test_label, B_test_label])
    return X_train, y_train, X_test, y_test

train_data, train_labels, test_data, _ = make_train_test_sets()

そして訓練データ train_data と train_labels で SVM の学習を行います。

from sklearn.svm import SVC

train_kernel = calculate_kernel(feature_map, train_data)

model = SVC(kernel='precomputed')
model.fit(train_kernel, train_labels);

これで学習できたはずなので、テストデータ test_data で検証を行います。

import matplotlib

test_kernel = calculate_kernel(feature_map, train_data, test_data)
pred = model.predict(test_kernel)

fig = matplotlib.pyplot.figure(figsize=(5, 5))
ax = fig.add_subplot(111)
ax.set_aspect('equal')
ax.set_title("Predicted data classification")
ax.set_ylim(-2, 2)
ax.set_xlim(-2, 2)
for (x, y), pred_label in zip(test_data, pred):
    c = 'C0' if pred_label == 0 else 'C3'
    ax.add_patch(matplotlib.patches.Circle((x, y), radius=.01,
                 fill=True, linestyle='solid', linewidth=4.0,
                 color=c))
matplotlib.pyplot.grid()
matplotlib.pyplot.show()

<Figure size 360x360 with 1 Axes>

大体それっぽい感じの結果になったのではないでしょうか。

え？これだけ？全然嬉しくないんだけど・・・

そうですね・・・。簡単な問題を解いたので何も嬉しくないですね。多項式カーネルや RBF カーネルで十分ですね。実際に自作カーネルを使うにしても素朴な特徴写像とグラム行列で十分です。これについてはご興味がありましたら拙記事 カーネル SVM を眺めてみる をごらんください。ただ、もっと複雑なデータセットの場合だと RBF カーネルで分類させると以下のようになることがあります。

あまり嬉しくない感じですね。

• こういう古典的に扱うのが難しいような複雑なデータセットが見つけられた場合に
• 古典計算で容易に近似できないような量子カーネルをうまく見つけられると

良い分類ができて 量子優位性が出るんじゃないかな？ というのが量子カーネル SVM の主張のようです。ということで機械的に適用すれば量子優位性が得られるような話ではなく、適したデータセットと適したカーネルによって量子優位性が得られることが期待できるかもしれない というお話となります。

まとめ

かなり無理やり詰め込んだので果たして当初の目的を達成できたか不安はありますが、「簡単なデータセットとフルスクラッチで書いた簡単な量子回路で問題を解いてみた」ということをしてみました。

但し、今回扱ったようなカーネルは初期の QSVM で提案されたものであり、最近の論文によると汎化性能が良くないようです。また初期のカーネルではデータサイズに対して計算量が 2 次的に増大するというスケーラビリティに関する問題もあるようです。これらに対応すべく新しく提案されたものについては例えば論文 [N-T-Y] などが参考になるかもしれません。ここで概要をまとめられると良かったのですが、すみません、ちょうど知ったばかりのやつでまだ読めていませんので詳細は割愛します。今記事の続きは論文 [N-T-Y] や文献 [QISKIT] で理解を深めていただければと思います。

参考文献

• [QISKIT] Quantum feature maps and kernels, Qiskit textbook
• [H-C-T-H-K-C-G] Vojtech Havlicek, Antonio D. Córcoles, Kristan Temme, Aram W. Harrow, Abhinav Kandala, Jerry M. Chow, Jay M. Gambetta, Supervised learning with quantum enhanced feature spaces, arXiv, 2018
• [N-T-Y] Kouhei Nakaji, Hiroyuki Tezuka, Naoki Yamamoto, Deterministic and random features for large-scale quantum kernel machine, arXiv, 2022
• [PML] S. Raschka, Python機械学習プログラミング, インプレス, 2020

拙記事

本記事を書く前に色々と書いていた記事がありますので参考までにリンクを掲載しておきます。本記事はこれらを思いっきり凝縮して単一の記事にしたものとも言えます。式や画像はこれら既存の記事から切り貼りしました。

• 最適化について考える (1) — Lagrange の未定乗数法
• 最適化について考える (2) — SVM
• カーネル SVM を眺めてみる
• 量子カーネル SVM で遊んでみる — Qiskit
• 量子カーネル SVM で遊んでみる — Blueqat