量子モンテ カルロ (QMC) シミュレーション

この例では、MATLAB® で量子モンテ カルロ (QMC) シミュレーションを使用して、確率変数の関数の平均を計算する方法を示します。オプション価格決定からマクロ経済のストレス テストまで、金融と経済の幅広いタスクで、モンテ カルロ シミュレーションが利用されています。この例では計算効率については検討しませんが、研究により、QMC は古典モンテ カルロ法に比べて二次的な高速化を実現することがわかっています。

この例では [1] に従って、確率変数の関数の平均を計算するために QMC を適用します。具体的には、基礎となる正規分布の三角関数の平均を計算します。この計算は、多くの実用的な用途に一般化されます。たとえば、確率分布が原資産の価格を表す場合、関数はその資産のオプションの価格にすることができます。

問題の定式化

正規分布に従う確率変数 "x" があり、評価する関数が以下であるとします。

$$\mathit{f}\left(\mathit{x}\right)={\sin }^2\left(\mathit{x}\right)$$

$\mathit{f}\left(\mathit{x}\right)$ の期待値は解析的に $\mu =\frac{\sinh \left(1\right)}{\mathit{e}}$ によって与えられます。

$\mu$ の値を計算します。

AnalyticMean = sinh(1)/exp(1)

AnalyticMean = 
0.4323

古典モンテカルロ法

比較のために、古典モンテカルロ法で平均 $\mu$ を計算します。その手順は以下のとおりです。

1. "x" の基となる分布から点のサンプルを生成します。

2. サンプル内の各点における関数値を計算します。

3. 関数値のサンプル平均を計算します。

func = @(x) sin(x).^2;
numSamples = 1000;
sampleData = randn(numSamples,1);
funcValues = func(sampleData);
MCMean = mean(funcValues)

MCMean = 
0.4210

サンプル、関数値、および計算された平均をプロットします。

tiledlayout(1,2)
nexttile
h1 = histogram(sampleData);
title("Sample Points")
nexttile
h2 = histogram(funcValues);
hold on
xline(MCMean,'--','MCMean');
hold off
title("Classic Monte Carlo","Mean of Function Values")

QMC 法

まず、量子シミュレーションのいくつかのパラメーターを定義します。確率分布レジスタの量子ビット数を 5 として定義し、推定レジスタの量子ビット数を 6 として定義します。次に、グリッドで確率分布および関数を離散化し、関数の離散平均を計算します。

m = 5;  % Number of probability qubits
n = 6;  % Number of estimation qubits
M = 2^m;
N = 2^n;
probQubits = 1:m; 

x_max = pi;
x = linspace(-x_max,x_max,M)';
p = normpdf(x);
p = p./sum(p);

DiscreteMean = sum(func(x).*p)

DiscreteMean = 
0.4326

[1] で説明されている QMC 法は、以下の手順で構成されています。

1. 確率分布を m 個の量子ビットの量子レジスタに読み込む。

2. 確率変数の値を値量子ビットに符号化する。

3. n 個の量子ビットのレジスタとともに振幅推定を使用して、値量子ビットの状態が $\left|1\right.anglebracketright$ である確率を推定します。

振幅推定には他の手法 (反復的な振幅推定など) もありますが、この例では量子位相推定を使用します。

以下の回路図に、この手法の概要を示します。F ブロックは確率分布を読み込み、確率変数を値量子ビットに符号化します。Q ブロックと逆量子フーリエ変換 (QFT) は量子位相推定に使用されます。

1. 確率分布の読み込み

確率分布を m 個の量子ビットに読み込みます。[1] では量子変分回路を使用して確率分布を近似していますが、この例では確率を回路に直接読み込みます。initGate 関数を使用して回路を初期化します。

A = initGate(probQubits,sqrt(p));
A.Name = 'Ainit';

レジスタをシミュレーションして、確率分布が正しく読み込まれたことを確認します。

circ = quantumCircuit(A);
state = simulate(circ);
figure
histogram(state,probQubits)
title("Simulated Probabilities")

2. 確率変数の符号化

離散化された確率変数関数を値量子ビットに符号化します。

valueQubit = m + 1;                        % Value qubit index
estQubits = (valueQubit + 1):(n + m + 1);  % Estimation qubits index

anglesR = 2*asin(sqrt(func(x)));

R = inv(ucryGate(probQubits,valueQubit,anglesR));
R.Name = "R";

値量子ビットの状態が $\left|1\right.anglebracketright$ である確率を計算します。F 回路の目的は、確率変数関数の期待値を評価することです。

F = compositeGate([A; R],[probQubits valueQubit],Name="F");
sF = simulate(quantumCircuit(F));
probability(sF,valueQubit,"1")

ans = 
0.4326

3. 量子位相推定

量子位相推定を使用して値量子ビットの値を推定します。

回路図の制御 Q ブロックは、[1] に概説されているように、複数の複合ゲートを用いて実装されています。この回路要素を表す複合ゲートを作成し、回路をプロットして含まれているゲートを確認します。

controlQubit = valueQubit+1; % Re-mapped to estQubits in the loop below

cV = czGate(controlQubit,valueQubit); % For simplicity, do not wrap in a named CompositeGate

cZgates = [xGate([probQubits valueQubit])
          hGate(probQubits(1))
          mcxGate([probQubits(2:end) valueQubit controlQubit],probQubits(1),[])
          hGate(probQubits(1))
          xGate([probQubits valueQubit])];

cZ = compositeGate(cZgates,1:controlQubit,Name="cZ");

% Put the gates together to form controlled Q:

cQ = compositeGate([cV; inv(F); cZ; F; cV; inv(F); cZ; F],...
                   [probQubits valueQubit controlQubit],Name='cQ');

QPE = [];
for ii=1:n
    % Repeat the cQ gate as needed and map to the current control qubit:
    cQrep = compositeGate(repmat(cQ,2^(n-ii),1),...
                          [probQubits valueQubit estQubits(ii)],...
                          Name="cQ"+2^(n-ii));

    QPE = [QPE; cQrep];
end

plot(cQ)

最終回路のシミュレーション

さまざまなコンポーネント回路を使用して最終回路を組み立て、結果をプロットします。名前と値の引数 QubitBlocks を指定して、2 つのレジスタおよび単一の値量子ビットを表示します。

circ = quantumCircuit([F;
                       hGate(estQubits); 
                       QPE; 
                       inv(qftGate(estQubits))]);
plot(circ,"QubitBlocks",[m 1 n])

回路をシミュレーションします。

state = simulate(circ);

結果の比較

位相 $\theta$ を用いて、解析平均は次のように与えられます。

$$\mu =\frac{1-\cos \left(\pi \theta \right)}{2}$$

$\mu$ の解析値を使用して $\theta$ について解きます。

AnalyticPhase = acos(-2*AnalyticMean + 1)/pi

AnalyticPhase = 
0.4568

解析位相と量子位相の推定を比較します。

[states,probabilities] = querystates(state,estQubits);
plot(bin2dec(states)/N,probabilities)
xlim([.35 .65])
xline(AnalyticPhase,'--',"AnalyticPhase")
xlabel("\theta")
ylabel("Probability")
title("Phase Estimation with QMC")

推定位相は振幅の最初の最大値です。推定位相を使用して関数の平均値を計算し、量子値を解析値、古典モンテカルロ法の値、および離散値と比較します。

y_ind = find(abs(probabilities - max(probabilities)) < 1e3*eps,1);
phase_estimated = bin2dec(states(y_ind))/N;
QMCMean = (1 - cos(pi * phase_estimated)) / 2;

ResultTable = table(AnalyticMean,MCMean,DiscreteMean,QMCMean)

ResultTable=1×4 table
    AnalyticMean    MCMean     DiscreteMean    QMCMean
    ____________    _______    ____________    _______

      0.43233       0.42103      0.43264       0.42663

参考文献

[1] Priazhkina, S., V. Skavysh, D. Guala, & T. Bromley."Quantum Monte Carlo for Economics: Stress Testing and Macroeconomic Deep Learning."Bank of Canada working paper (2022). https://www.bankofcanada.ca/wp-content/uploads/2022/06/swp2022-29.pdf

[2] Brassard, G., P. Hoyer, M. Mosca, & A. Tapp."Quantum Amplitude Amplification and Estimation."Contemporary Mathematics 305 (2002): 53-74. https://arxiv.org/pdf/quant-ph/0005055

[3] Rebentrost, P., B. Gupt, & T. R. Bromley."Quantum Computational Finance: Monte Carlo Pricing of Financial Derivatives."Physical Review A, 98(2), 022321 (2018). https://arxiv.org/pdf/1805.00109

[4] Stamatopoulos, N., D. J. Egger, Y. Sun, C. Zoufal, R. Iten, N. Shen, & S. Woerner."Option Pricing Using Quantum Computers."Quantum 4, 291 (2020). https://arxiv.org/pdf/1905.02666

[5] Woerner, S., & D. J. Egger."Quantum Risk Analysis." npj Quantum Information 5(1), 15 (2019). https://arxiv.org/pdf/1806.06893

参考

quantumCircuit

トピック

• Introduction to Quantum Computing
• ローカル量子状態シミュレーション
• Graph Coloring with Grover's Algorithm
• Ground-State Protein Folding Using Variational Quantum Eigensolver (VQE)
• Solve XOR Problem Using Quantum Neural Network (QNN)

外部の Web サイト

• Quantum Monte Carlo with MATLAB