アニーリングを用いたクラスタリング#

このチュートリアルでは、アニーリングの応用の一例としてPyQUBOとOpenjijを利用したクラスタリングを取り上げます。

クラスタリング#

クラスタリングとは与えられたデータをnn個のクラスターに分けるというタスクです(nnは外部から与えられているとします)。簡単のため、今回はクラスター数が2を考えましょう。

クラスタリングのハミルトニアン#

今回は、以下のハミルトニアンを最小化することでクラスタリングを行います。

H=−∑i,j12di,j(1−σiσj) H = - \sum_{i, j} \frac{1}{2}d_{i,j} (1 - \sigma _i \sigma_j)

i,ji, jはサンプルの番号、di,jd_{i,j}は２つのサンプル間の距離、σi={−1,1}\sigma_i=\{-1,1\}は２つのクラスターのどちらかに属しているかを表すスピン変数です。
このハミルトニアンの和の各項は

  • σi=σj\sigma_i = \sigma_j のとき、0

  • σi≠σj\sigma_i \neq \sigma_j のとき、di,jd_{i,j}

となります。右辺のマイナスに注意すると、ハミルトニアン全体では「異なるクラスに属しているサンプル同士の距離を最大にする{σ1,σ2…}\{\sigma _1, \sigma _2 \ldots \}の組を選びなさい」という問題に帰着することがわかります。

必要なライブラリのインポート#

!pip install jijmodeling-transpiler

Requirement already satisfied: jijmodeling-transpiler in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (0.3.8)
Requirement already satisfied: numpy<1.25.0,>=1.17.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (1.24.2)
Requirement already satisfied: pyqubo<1.4.0,>=1.3.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (1.3.1)
Requirement already satisfied: jijmodeling>=0.10.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (0.10.19)
Requirement already satisfied: mip in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (1.15.0)
Requirement already satisfied: typeguard in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (2.13.3)
Requirement already satisfied: pydantic in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-transpiler) (1.10.6)
Requirement already satisfied: python-ulid<1.2.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling>=0.10.0->jijmodeling-transpiler) (1.1.0)
Requirement already satisfied: pandas<1.6.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling>=0.10.0->jijmodeling-transpiler) (1.5.3)
Requirement already satisfied: jijmodeling-schema==0.0.17 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling>=0.10.0->jijmodeling-transpiler) (0.0.17)
Requirement already satisfied: orjson<3.9.0,>=3.8.7 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling>=0.10.0->jijmodeling-transpiler) (3.8.8)
Requirement already satisfied: betterproto[compiler]==1.2.5 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (1.2.5)
Requirement already satisfied: grpclib in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (0.4.3)
Requirement already satisfied: stringcase in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (1.2.0)
Requirement already satisfied: jinja2 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (3.1.2)
Requirement already satisfied: black in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (23.1.0)
Requirement already satisfied: protobuf in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (4.22.1)
Requirement already satisfied: typing-extensions>=4.2.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pydantic->jijmodeling-transpiler) (4.5.0)
Requirement already satisfied: dwave-neal>=0.5.7 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (0.6.0)
Requirement already satisfied: dimod<0.13,>=0.9.14 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (0.12.4)
Requirement already satisfied: six>=1.15.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (1.16.0)
Requirement already satisfied: Deprecated>=1.2.12 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (1.2.13)
Requirement already satisfied: cffi==1.15.* in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from mip->jijmodeling-transpiler) (1.15.1)
Requirement already satisfied: pycparser in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from cffi==1.15.*->mip->jijmodeling-transpiler) (2.21)
Requirement already satisfied: wrapt<2,>=1.10 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from Deprecated>=1.2.12->pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (1.15.0)
Requirement already satisfied: dwave-samplers<2.0.0,>=1.0.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from dwave-neal>=0.5.7->pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (1.0.0)
Requirement already satisfied: pytz>=2020.1 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pandas<1.6.0->jijmodeling>=0.10.0->jijmodeling-transpiler) (2022.7.1)
Requirement already satisfied: python-dateutil>=2.8.1 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from pandas<1.6.0->jijmodeling>=0.10.0->jijmodeling-transpiler) (2.8.2)
Requirement already satisfied: networkx>=2.4.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from dwave-samplers<2.0.0,>=1.0.0->dwave-neal>=0.5.7->pyqubo<1.4.0,>=1.3.0->jijmodeling-transpiler) (2.8.8)
Requirement already satisfied: mypy-extensions>=0.4.3 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (1.0.0)
Requirement already satisfied: platformdirs>=2 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (3.1.1)
Requirement already satisfied: click>=8.0.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (8.1.3)
Requirement already satisfied: packaging>=22.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (23.0)
Requirement already satisfied: tomli>=1.1.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (2.0.1)
Requirement already satisfied: pathspec>=0.9.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from black->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (0.11.1)
Requirement already satisfied: multidict in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from grpclib->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (6.0.4)
Requirement already satisfied: h2<5,>=3.1.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from grpclib->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (4.1.0)
Requirement already satisfied: MarkupSafe>=2.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from jinja2->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (2.1.2)
Requirement already satisfied: hyperframe<7,>=6.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from h2<5,>=3.1.0->grpclib->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (6.0.1)
Requirement already satisfied: hpack<5,>=4.0 in /opt/hostedtoolcache/Python/3.9.16/x64/lib/python3.9/site-packages (from h2<5,>=3.1.0->grpclib->betterproto[compiler]==1.2.5->jijmodeling-schema==0.0.17->jijmodeling>=0.10.0->jijmodeling-transpiler) (4.0.0)

# ライブラリのインポート
import numpy as np
from matplotlib import pyplot as plt
import pandas as pd
from scipy.spatial import distance_matrix

import openjij as oj
import jijmodeling as jm
from jijmodeling.transpiler.pyqubo.to_pyqubo import to_pyqubo

JijModelingとOpenJijを用いたクラスタリング#

最初に、JijModelingを用いた上記のハミルトニアンを定式化します。JijModelingでは、スピン変数σi\sigma_iを扱えないので、バイナリ変数xix_iで書けるように、σi=2xi−1\sigma_i = 2x_i - 1という関係を用いて書き直します。

problem = jm.Problem("clustering")
d = jm.Placeholder("d", dim=2)
N = d.shape[0].set_latex('N')
x = jm.Binary("x", shape=(N))
i = jm.Element("i", (0, N))
j = jm.Element("j", (0, N))
problem += (
    -1 / 2 * jm.Sum([i, j], d[i, j] * (1 - (2 * x[i] - 1) * (2 * x[j] - 1)))
)
problem
Problem: clusteringmin⁡(−0.5)⋅∑i=0N−1∑j=0N−1di,j⋅(1−(2⋅xi−1)⋅(2⋅xj−1))xi0∈{0,1}\begin{alignat*}{4}\text{Problem} & \text{: clustering} \\\min & \quad \left( -0.5 \right) \cdot \sum_{ i = 0 }^{ N - 1 } \sum_{ j = 0 }^{ N - 1 } d_{i,j} \cdot \left( 1 - \left( 2 \cdot x_{i} - 1 \right) \cdot \left( 2 \cdot x_{j} - 1 \right) \right) \\& x_{i_{0}} \in \{0, 1\}\end{alignat*}

人工データの生成#

今回は人工的に二次元平面上の明らかに線形分離可能なデータを生成しましょう。

data = []
label = []
for i in range(100):
    # [0, 1]の乱数を生成
    p = np.random.uniform(0, 1)
    # ある条件を満たすときをクラス1、満たさないときを-1
    cls =1 if p>0.5 else -1
    # ある座標を中心とする正規分布に従う乱数を作成
    data.append(np.random.normal(0, 0.5, 2) + np.array([cls, cls]))
    label.append(cls)
# DataFrameとして整形
df1 = pd.DataFrame(data, columns=["x", "y"], index=range(len(data)))
df1["label"] = label

# データセットの可視化
df1.plot(kind='scatter', x="x", y="y")

<Axes: xlabel='x', ylabel='y'>
../../../_images/ed84cd418f8adf5fc1215262cbb8eaffa7d2641c2910a4fa15e4c373933b9b75.png 
instance_data = {'d':distance_matrix(df1, df1)}

OpenJijを使ってクラスタリング問題を解く#

数理モデルとデータができたので、早速、Openjijを使って問題を解いていきましょう。

pyq_obj, pyq_cache = to_pyqubo(problem, instance_data, {})
qubo, constant = pyq_obj.compile().to_qubo()
sampler = oj.SASampler()
response = sampler.sample_qubo(qubo)
result = pyq_cache.decode(response)

可視化をしてみましょう。

for idx in range(0, len(instance_data['d'])):
    if idx in result.record.solution["x"][0][0][0]:
        plt.scatter(df1.loc[idx]["x"], df1.loc[idx]["y"], color="b")
    else:
        plt.scatter(df1.loc[idx]["x"], df1.loc[idx]["y"], color="r")
plt.xlabel('x')
plt.ylabel('y')

Text(0, 0.5, 'y')
../../../_images/eaff20dc0ce5a7abf179c580ee503aa64bf2061b6e9e277c9d858396d3ee3179.png

赤と青のクラスに分類できていることがわかります。