scikit-qulacs チュートリアル#

scikit-qulacs とは#

量子回路シミュレータ Qulacs-Osaka (Qulacs の大阪大学バージョン) をバックエンドとした量子機械学習アルゴリズムのシミュレータです。scikit-learn の名前をもじっています。

インストール#

PyPI からインストールできます。

pip install skqulacs

※Python 3.8以降を使用してください。Python 3.7系の場合、PyPIインストールで旧版0.2.0がインストールされます。最新コードで実装されているチュートリアルがエラーになる場合があります。

Quantum neural network#

入力データ \({\boldsymbol x}\) に対して、学習パラメータ \({\boldsymbol \theta}\) を使って

\[y_{{\boldsymbol \theta}}({\boldsymbol x}) = \langle 0|U^\dagger({\boldsymbol \theta}, {\boldsymbol x}) O U({\boldsymbol \theta}, {\boldsymbol x})|0\rangle\]

の形でモデル \(y_{{\boldsymbol \theta}}({\boldsymbol x})\) を構築する手法を量子ニューラルネットワークと呼びます。scikit-qulacs では skqulacs.qnn モジュールに実装されています。分類問題をとく skqulacs.qnn.classifier と回帰問題を解く skqulacs.qnn.regressor の２つがあります。

`skqulacs.qnn.QNNClassifier` の使い方#

例として、有名な iris データセットを分類する量子ニューラルネットワークを組んでみます。まず、以下のコードによって scikit-learn から iris データセットを取り出し、訓練データとテストデータに分けます。

import pandas as pd
from sklearn import datasets
from sklearn.model_selection import train_test_split

iris = datasets.load_iris()
df = pd.DataFrame(iris.data, columns=iris.feature_names)
x = df.loc[:, ["petal length (cm)", "petal width (cm)"]]
x_train, x_test, y_train, y_test = train_test_split(
    x, iris.target, test_size=0.25, random_state=0
)
x_train = x_train.to_numpy()
x_test = x_test.to_numpy()

変数はそれぞれ

x_train: 訓練データの入力
y_train: 訓練データのラベル
x_test: テストデータの入力
y_test: テストデータのラベル

です。

次に簡単な量子回路によってこれらの分類を行ってみます。

from skqulacs.circuit.pre_defined import create_qcl_ansatz
from skqulacs.qnn import QNNClassifier
from skqulacs.qnn.solver import Bfgs

nqubit = 5                      # qubitの数。必要とする入力、出力の次元数以上が必要である。
c_depth = 3                     # circuitの深さ。ニューラルネットワークでの隠れ層に対応する。
time_step = 1.                  # qcl_ansatz に必要となるパラメータ。横磁場イジングモデルによって発展させる時間を決める。
num_class = 3                   # 分類数（ここでは3つの品種に分類）
solver = Bfgs()                   # アルゴリズム。ほかには、`Adam()` や `NelderMead()` が使える。
maxiter = 20                    # ループの最大。これが多いほど、正確になるが、時間がかかる。
circuit = create_qcl_ansatz(nqubit, c_depth, time_step)     # LearningCircuitを作る
qcl = QNNClassifier(circuit, num_class, solver)     # モデル構築を行う
opt_loss, opt_params = qcl.fit(x_train, y_train, maxiter)                          # 学習
print("trained parameters", opt_params)
print("loss", opt_loss)

trained parameters [ 3.3077862   2.20610834  0.18317852  2.57427866  2.1731261   5.16437333
90379975  4.92301026  3.78621276  5.01264629  5.39612169  5.06373429
36012806  4.08739844  5.47151574  2.70898056  3.1735183   2.66824882
2683309   4.5094121   6.64751394 -0.05754961  6.58882916  4.47236961
6969534   6.44925444  0.82493735  5.89706655  0.36732051  4.13386748
63948049  2.33409167  4.17240084  0.31197252  2.18623163  0.46591338
45247138  4.11809259  0.99494374  4.61047697  6.022204    3.2083882
44504131  2.8614801   5.06963927]
loss 0.3253091864984035

テストデータと比べて結果を出力します。

from sklearn.metrics import f1_score
y_pred = qcl.predict(x_test)
print(f1_score(y_test, y_pred, average="weighted"))

0.9473684210526315

うまく分類できていることがわかります。

`skqulacs.qnn.QNNRegressor` の使い方#

回帰問題を解くには QNNRegressor を使います。ここでは \(y=\sin(\pi x)\) を学習させてみます。まずデータセットを作ります。

import numpy as np
from numpy.random import default_rng

def generate_noisy_sine(x_min, x_max, num_x):
    rng = default_rng(0)
    x_train = [[rng.uniform(x_min, x_max)] for _ in range(num_x)]
    y_train = [np.sin(np.pi*x[0]) for x in x_train]
    mag_noise = 0.01
    y_train += mag_noise * rng.random(num_x)
    return np.array(x_train), np.array(y_train)
  
x_min = -1.0
x_max = 1.0
num_x = 80
x_train, y_train = generate_noisy_sine(x_min, x_max, num_x)
x_test, y_test = generate_noisy_sine(x_min, x_max, num_x)

作成したデータに対して学習を行います。

from skqulacs.qnn import QNNRegressor

n_qubit = 4
depth = 6
time_step = 0.5
solver = Bfgs()
maxiter= 30
circuit = create_qcl_ansatz(n_qubit, depth, time_step, 0)
qnn = QNNRegressor(circuit, solver)
opt_loss, opt_params = qnn.fit(x_train, y_train, maxiter)
print("trained parameters", opt_params)
print("loss", opt_loss)

trained parameters [ 4.25798056  1.92999783  0.31115392  0.49212499  4.80907853  6.06343962
00173116  4.9030698   3.43327385  5.81135678  5.07675043  0.49105851
42057489 -0.14168464  4.58566892  1.36858878  5.36424976  3.55819331
73062387  2.57575324  0.13021265  1.13880549  4.24069209  3.98487665
57694613  2.3916734   6.35641985  6.29467416  4.37029016  4.04927805
49274116  2.73557408  0.63853703  4.45377244  3.03068676  2.05579799
04866402  5.44965241  5.85460933  2.23874886  3.3939702   2.04538289
79038304  1.9171457   2.36846662  5.80016494  1.49691542  4.21257824
37117935  5.12489574  4.95969291  1.41616941  5.50035696  0.29056759
13739014  0.89954726  2.93673013  5.33719954  1.33472466  0.45390472
59255057  1.20163753  0.38793178  3.64633592  1.87676315  4.22226834
25359251  5.91947125  2.29405484  0.66284639  3.95280309  5.82548387]
loss 0.001734007213230344

訓練したモデルとテストデータを並べてプロットしてみます。

import matplotlib.pyplot as plt

y_pred = qnn.predict(x_test)
plt.plot(x_test, y_test, "o", label="Test")
plt.plot(np.sort(np.array(x_test).flatten()), np.array(y_pred)[np.argsort(np.array(x_test).flatten())], label="Prediction")
plt.legend()
plt.show()

../_images/46b5b48546faa84faaf45862bc9e3690f3ca4a5d2c4d2616d5207210ef04c401.png

よく学習できていることがわかります。

`skqulacs.qnn.generation` の使い方#

生成モデルを作ります。生成モデルとは、教師なし学習の一種で、入力データの分布を学習し、それに応じた新しいサンプルを生成するモデルです。

ここでは、例として、ガウス関数の和のモデルを学習させてみます。

from math import exp

import numpy as np
import pytest

from skqulacs.circuit import create_farhi_neven_ansatz
from skqulacs.qnn import QNNGeneretor
from skqulacs.qnn.solver import Bfgs
n_qubit = 5
depth = 10
circuit = create_farhi_neven_ansatz(n_qubit, depth)

qnn = QNNGeneretor(circuit, Bfgs(), "same", 4.0, 5)

#prob_listつまり、確率分布を生成する
prob_list = np.zeros(32)
ua = 32 * 2 / 7
ub = 32 * 5 / 7
v = 32 * 1 / 8
prob_sum = 0
for i in range(32):
    prob_list[i] = exp(-(ua - i) * (ua - i) / (2 * v * v)) + exp(
        -(ub - i) * (ub - i) / (2 * v * v)
    )
    prob_sum += prob_list[i]

for i in range(32):
    prob_list[i] /= prob_sum

#prob_list が、理想的な確率分布になる。

#prob_listの確率分布をもとに、10000個のデータを作る
datas = np.random.choice(a=range(32), size=10000, p=prob_list)


maxiter = 120
qnn.fit(datas, maxiter)

data_param = qnn.predict()

/home/watle/pros/scikit-qulacs/.venv/lib/python3.8/site-packages/scipy/optimize/linesearch.py:167: ComplexWarning: Casting complex values to real discards the imaginary part
  stp, phi1, derphi1, task = minpack2.dcsrch(alpha1, phi1, derphi1,

plt.plot(list(range(0,32)), data_param)
plt.plot(list(range(0,32)), prob_list)

/home/watle/pros/scikit-qulacs/.venv/lib/python3.8/site-packages/matplotlib/cbook/__init__.py:1298: ComplexWarning: Casting complex values to real discards the imaginary part
  return np.asarray(x, float)

[<matplotlib.lines.Line2D at 0x7fec5fd9fd60>]

../_images/c9df8c4e05899daf09da435523cac1056e1d4dfefc75c27ee996fc2729d4ff4f.png

オレンジが学習すべきものの頻度の折れ線グラフ、青が学習した頻度の折れ線グラフです。

うまく学習できています。

量子カーネル法#

入力データ \({\boldsymbol x}\) を量子状態 \(|\psi({\boldsymbol x})\rangle\) にマップし、この量子状態間の内積で定義されるカーネル関数

\[K({\boldsymbol x},{\boldsymbol x}')=|\langle \psi({\boldsymbol x})| \psi({\boldsymbol x}')\rangle|^2\]

を使ったカーネル法を量子カーネル法と呼びます。scikit-qulacs には、support vector classifier skqulacs.qsvm.QSVC と support vector regressor skqulacs.qsvm.QSVR の2つが実装されています。

`skqulacs.qsvm.QSVC` の使い方#

ここでは IBM のグループが論文 Supervised learning with quantum-enhanced feature spaces で提案した \(|\psi({\boldsymbol x})\rangle\) を使った分類を行います。この状態を作る量子回路は skqulacs.circuit.create_ibm_embedding_circuit によって作り出せます。

まずデータセットを準備します。

iris = datasets.load_iris()
df = pd.DataFrame(iris.data, columns=iris.feature_names)
x = df.loc[:, ["petal length (cm)", "petal width (cm)"]]
x_train, x_test, y_train, y_test = train_test_split(
    x, iris.target, test_size=0.25, random_state=0
)
x_train = x_train.to_numpy()
x_test = x_test.to_numpy()

回路を作成し、QSVC を訓練、予測値を出力させます。

from skqulacs.circuit import create_ibm_embedding_circuit
from skqulacs.qsvm import QSVC
n_qubit = 4                         # x_train の次元数以上必要。あまり小さいと結果が悪くなる。
circuit = create_ibm_embedding_circuit(n_qubit)
qsvm = QSVC(circuit)
qsvm.fit(x_train, y_train)
y_pred = qsvm.predict(x_test)

評価してみます。

print(f1_score(y_test, y_pred, average="weighted"))

0.8978070175438595

qnn よりも悪い結果となりました。回路形がこのデータセットに適していないからだと考えられます。

`skqulacs.qsvm.QSVR` の使い方#

同様の回路を使い、回帰を行ってみます。QNN の場合と同様に、\(\sin(\pi x)\) を学習します。まずデータセットを作ります。

x_min = -1.0
x_max = 1.0
num_x = 1000
x_train, y_train = generate_noisy_sine(x_min, x_max, num_x)
x_test, y_test = generate_noisy_sine(x_min, x_max, num_x)

学習します。

from skqulacs.qsvm import QSVR
n_qubit = 4                         # x_train の次元数以上必要。あまり小さいと結果が悪くなる。
circuit = create_ibm_embedding_circuit(n_qubit)
qsvm = QSVR(circuit)
qsvm.fit(x_train, y_train)
y_pred = qsvm.predict(x_test)

学習結果をプロットしてみます。

plt.plot(x_test, y_test, "o", label="Test")
plt.plot(np.sort(np.array(x_test).flatten()), np.array(y_pred)[np.argsort(np.array(x_test).flatten())], label="Prediction")
plt.legend()
plt.show()

../_images/880ffe114c99d4dd4e5a10f6391db33d20912414da61c5b2dcf47667f01ed7b2.png

うまく学習できています。

用意されている量子回路#

scikit-qulacs では代表的な論文で発表されている機械学習用量子回路を手軽に使えるようにするため、量子回路の実装に努めています。機械学習用量子回路は LearningCircuit クラスのインスタンスとして実装されていて、skqulacs.circuit モジュールにある create_...() という関数を呼ぶことにより作成可能です。

代表的なansatzとして、

create_farhi_neven_ansatz(n_qubit: int, c_depth: int, seed: Optional[int] = None): arXiv:1802.06002
create_ibm_embedding_circuit(n_qubit: int): arXiv:1804.11326

が用意されています。

他に利用可能な量子回路のリストについてはこのページをご覧ください。

基本的には QNNClassifier などの引数として与えることを推奨しますが、自前で新たなアルゴリズムを作りたいときなど、これらの量子回路を直接扱いたいときもあると思います。以下に使用例を掲載します。

以下のコードで QCL ansatz を作成します。

n_qubits = 4
c_depth = 4
time_step = 1.
ansatz = create_qcl_ansatz(n_qubits, c_depth, time_step)

LearningCircuit にセットされている訓練パラメータは get_parameters() によって取り出せます。

parameters = ansatz.get_parameters()
print(len(parameters), parameters)      # show parameters in the circuit. Parameters are chosen randomly within create_qcl_ansatz()

48 [2.9646519753472966, 2.4054163630603496, 5.8984153363073375, 2.512801032804015, 5.964913519260788, 3.8936080016830195, 5.923572589632018, 4.9481732525951205, 1.568563629432675, 5.554018967676112, 0.028389796241806445, 2.897527605841249, 5.267355240624224, 3.735313799160573, 2.2661746616128897, 3.5452804132198086, 0.43246053635497056, 2.5730108963703193, 1.974663878510053, 1.9654655970918322, 3.967587987184092, 0.8120259125764905, 6.047592885353966, 0.5745413307337336, 0.9710274726615847, 0.4695106940105736, 2.341674095381348, 0.802699963876766, 0.8310916316734924, 5.012384620423008, 0.018233484756863608, 4.2472370241118345, 5.1721310224316355, 5.021016841850289, 3.1032612340282735, 6.090572550277523, 4.484895234385721, 2.0890853870812838, 1.8127693347262788, 3.7158089135425, 3.7028938225792483, 3.3726572966673407, 2.699611657167366, 1.6621937249934382, 1.1061035348895243, 1.5957931515488104, 0.26282289161037825, 2.4484720212561286]

update_parameters() を使うと、パラメータを更新できます。以下では全てのパラメータをゼロにセットします。

ansatz.update_parameters(np.zeros(len(parameters)))
print(ansatz.get_parameters())

[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

run(x) を使うと、入力データ \({\boldsymbol x}\) に対して、初期状態 \(|0\rangle\) に量子回路を作用させたときの出力状態が qulacs.QuantumState として返されます。QCL ansatz は n_qubits 次元までの入力を扱うことのできる ansatz なので、今回は 2 次元入力を与えます。

多くの回路では、データ次元数をaとして、それがn_qubits以下の場合、回路のx番目には、n_qubits % a 番目の入力が入ります。

state = ansatz.run([np.pi/2, np.pi/2])
print(state.get_vector())

[ 0.30869612+0.13283004j -0.0759879 -0.06176279j -0.21316566-0.1438611j
21877642+0.06614673j -0.33638863+0.12504632j  0.24610584+0.23178805j
05948781+0.05711437j  0.03556247-0.08914448j  0.03556247-0.08914448j
05948781+0.05711437j  0.24610584+0.23178805j -0.33638863+0.12504632j
21877642+0.06614673j -0.21316566-0.1438611j  -0.0759879 -0.06176279j
30869612+0.13283004j]

ユーザー定義の量子回路を学習に利用する#

ここでは自分でオリジナルの量子回路を定義し、それを学習に利用する方法を説明します。scikit-qulacs では、LearningCircuit クラスによって量子回路を表現しています。これは qulacs の ParametricQuantumCircuit のラッパーとなっていて、入力データと訓練される回路パラメータを一元管理し、それらをすべて qulacs の parametric gate で扱えるようにしたものです。

LearningCircuit クラスは qulacs の ParametricQuantumCircuit のうちの主要な add_.._gate() 系のメソッドを機能そのままにラップしています。それに加えて、ParametricQuantumCircuit を拡張する機能として、以下のメソッドが追加されています。

add_input_RX_gate(index, input_func), add_input_RY_gate(index, input_func), add_input_RZ_gate(index, input_func): 入力データを入れるためのゲートを追加します。引数にはそれぞれ以下のものを与えます。
- index: そのゲートを作用させる量子ビット
- input_func: 入力データ \({\boldsymbol x}\) -> 回転ゲートの回転角 \(\phi = \phi({\boldsymbol x})\) の変換を行う関数。回路実行時には input_func(x) が角度にセットされます。
add_parametric_input_RX_gate(index, parameter, input_func) (RY, RZ も同様) :入力データと訓練パラメータが同時に含まれているようなゲートを追加します。
- index: そのゲートを作用させる量子ビット
- parameter: 実数の訓練パラメータ \(\theta\)
- input_func: 入力データ \({\boldsymbol x}\), 訓練パラメータ \(\theta\) -> 回転ゲートの回転角 \(\phi = \phi(\theta, {\boldsymbol x})\) の変換を行う関数。回路実行時には input_func(theta, x) が角度にセットされます。

例えば以下のコードは、各量子ビットに \({\boldsymbol x} = \{x_i\}_{i=0}^3\) の各成分を \(R_y\) ゲートで入力したあと、各量子ビットを訓練可能な \(R_y\) ゲートで変換し、最後に全ての量子ビットにアダマールゲートを一度だけ作用させる LearningCircuit を構成するコードです。(例のためだけに書いていて、特に意味は無い回路です。)

from skqulacs.circuit import LearningCircuit
n_qubit = 4
circuit = LearningCircuit(n_qubit)
for i in range(n_qubit):
    circuit.add_input_RY_gate(i, lambda x: x[i%2])  # input 2 dimensional x into the circuit by RY rotation
for i in range(n_qubit):
    circuit.add_parametric_RY_gate(i, 0.)           # first argument is qubit index, second argument is parameter
for i in range(n_qubit):
    circuit.add_H_gate(i)

run() メソッドによって初期状態 \(|0\rangle\) にこの回路を作用させたときの量子状態が返されます。このメソッドは入力データ \({\boldsymbol x}\) を引数として取ります。上の回路は 2 次元の \({\boldsymbol x}\) を取れるように設計したので、2次元の入力を与えます。

x = [np.pi/2,np.pi/2]
state = circuit.run(x)
print(state.get_vector())

[-6.93889390e-18+0.j -1.96261557e-17+0.j -2.08166817e-17+0.j
 -2.77555756e-17+0.j -2.08166817e-17+0.j  1.96261557e-17+0.j
 -6.93889390e-18+0.j  0.00000000e+00+0.j -2.08166817e-17+0.j
  1.96261557e-17+0.j -6.93889390e-18+0.j  2.77555756e-17+0.j
 -6.93889390e-18+0.j  5.88784672e-17+0.j  9.02056208e-17+0.j
  1.00000000e+00+0.j]

予想どおりの出力が得られることが確認できます。

scikit-qulacs チュートリアル

Contents