참고

이 페이지는 tutorials/finance/09_credit_risk_analysis.ipynb 에서 생성되었다.

IBM 퀀텀 랩 `IBM Quantum lab 에서 대화식으로 실행하시오.

신용 위험 분석

소개

This tutorial shows how quantum algorithms can be used for credit risk analysis. More precisely, how Quantum Amplitude Estimation (QAE) can be used to estimate risk measures with a quadratic speed-up over classical Monte Carlo simulation. The tutorial is based on the following papers:

• Quantum Risk Analysis. Stefan Woerner, Daniel J. Egger. [Woerner2019]

• Credit Risk Analysis using Quantum Computers. Egger et al. (2019) [Egger2019]

A general introduction to QAE can be found in the following paper:

• Quantum Amplitude Amplification and Estimation. Gilles Brassard et al.

The structure of the tutorial is as follows:

1. Problem Definition

2. Uncertainty Model

3. Expected Loss

4. Cumulative Distribution Function

5. Value at Risk

6. Conditional Value at Risk

import numpy as np
import matplotlib.pyplot as plt

from qiskit import QuantumRegister, QuantumCircuit, Aer, execute
from qiskit.circuit.library import IntegerComparator
from qiskit.aqua.algorithms import IterativeAmplitudeEstimation

문제점 정의

이 튜토리얼에서는 :math:의 K 자산 포트폴리오의 신용 위험을 분석하고자 합니다. 모든 자산 ( :math:’k’) 의 디폴트 확률은 *Gaussian 조건부 독립 * 모델을 따르며, 즉 표준 정규 분포에 이어서 잠재적인 무작위 변수 ( :math:’Z’) 로부터 샘플링된 값 :math:’z’ 가 제공되며, 자산 :math:’k’ 의 디폴트 확률은 다음과 같이 주어진다.

\[(z) = F\left (\frac {F={-1}(p_k ^ 0)-\sqrt {\rho_k} z} {\sqrt {1- \rho_k}} \right)\]

여기서 \(F\)rho_k`는 :math:`Z`에 대한 자산 :math:`k`의 기본 확률에 대한 민감도이다. 따라서 :math:`Z`의 구체적인 실현이 주어지면 개별 기본 이벤트는 서로 독립적 인 것으로 간주된다.

우리는 전체 손실의 위험 측정을 분석하는데 관심이 있습니다.

\[L = \sum_ {k= 1} ^ K \lambda_k X_k (Z)\]

where \(\lambda_k\) denotes the loss given default of asset \(k\), and given \(Z\), \(X_k(Z)\) denotes a Bernoulli variable representing the default event of asset \(k\). More precisely, we are interested in the expected value \(\mathbb{E}[L]\), the Value at Risk (VaR) of \(L\) and the Conditional Value at Risk of \(L\) (also called Expected Shortfall). Where VaR and CVaR are defined as

\[\text{VaR}_{\alpha}(L) = \inf \{ x \mid \mathbb{P}[L <= x] \geq 1 - \alpha \}\]

with confidence level \(\alpha \in [0, 1]\), and

\[\text{CVaR}_ {\alpha} (L) = \mathbb{E}[L \mid L \geq \text{VaR}_ {\alpha} (L)].\]

고려된 모델에 대한 자세한 내용은 예를 들어 신용 위험에 대한 규제 자본 모델링을 참조하십시오. 마렉 루트코스키 (Marek Rutkowski), 실비오 타카 (Silvio Tarca)

문제는 다음 매개 변수에 의해 정의된다.-\(Z\), 즉, Z는 \(2^{n_z}\){-z_{max}, …, + z_{max}}`-각 자산의 기본 기본 확률 \(p_0^k \in (0, 1)\), \(k = 1, ..., K\)-다음에 대한 기본 확률의 민감도 \(Z\), \(\rho_k \in [0, 1)\)-자산에 대한 기본 손실 \(k\), \(\lambda_k\)alpha in [0, 1]`.

# set problem parameters
n_z = 2
z_max = 2
z_values = np.linspace(-z_max, z_max, 2**n_z)
p_zeros = [0.15, 0.25]
rhos = [0.1, 0.05]
lgd = [1, 2]
K = len(p_zeros)
alpha = 0.05

불확실성 모델

이제 불확실성 모델을로드하는 회로를 구성합니다. 이것은 표준 정규 분포에 따라 : math :Z`를 나타내는 : math :`n_z 큐 비트의 레지스터에 양자 상태를 생성함으로써 달성 할 수 있다. 그런 다음이 상태는 \(K\) 큐 비트의 두 번째 큐 비트 레지스터에서 단일 큐 비트 Y 회전을 제어하는 ​​데 사용된다. 여기서 큐 비트 \(k\) 상태는 기본 이벤트를 나타낸다. 자산의 \(k\). 결과 양자 상태는 다음과 같이 쓸 수 있다.

\[ |\Psi\rangle = \sum_{i=0}^{2^{n_z}-1} \sqrt{p_z^i} |z_i \rangle \bigotimes_{k=1}^K \left( \sqrt{1 - p_k(z_i)}|0\rangle + \sqrt{p_k(z_i)}|1\rangle\right),\]

where we denote by \(z_i\) the \(i\)-th value of the discretized and truncated \(Z\) [Egger2019].

from qiskit.finance.applications import GaussianConditionalIndependenceModel as GCI
u = GCI(n_z, z_max, p_zeros, rhos)

u.draw()

     ┌───────┐┌─────────┐┌─────────┐
q_0: ┤0      ├┤0        ├┤0        ├
     │  P(X) ││         ││         │
q_1: ┤1      ├┤1 LinRot ├┤1        ├
     └───────┘│         ││  LinRot │
q_2: ─────────┤2        ├┤         ├
              └─────────┘│         │
q_3: ────────────────────┤2        ├
                         └─────────┘

이제 시뮬레이터를 사용하여 \(| \Psi \rangle\)mathbb{E} [L]`-PDF 및 CDF의 해당 정확한 값을 계산한다. \(L\)-위험 값 \(VaR (L)\)

# run the circuit and analyze the results
job = execute(u, backend=Aer.get_backend('statevector_simulator'))

# analyze uncertainty circuit and determine exact solutions
p_z = np.zeros(2**n_z)
p_default = np.zeros(K)
values = []
probabilities = []
num_qubits = u.num_qubits
for i, a in enumerate(job.result().get_statevector()):

    # get binary representation
    b = ('{0:0%sb}' % num_qubits).format(i)
    prob = np.abs(a)**2

    # extract value of Z and corresponding probability
    i_normal = int(b[-n_z:], 2)
    p_z[i_normal] += prob

    # determine overall default probability for k
    loss = 0
    for k in range(K):
        if b[K - k - 1] == '1':
            p_default[k] += prob
            loss += lgd[k]
    values += [loss]
    probabilities += [prob]

values = np.array(values)
probabilities = np.array(probabilities)

expected_loss = np.dot(values, probabilities)

losses = np.sort(np.unique(values))
pdf = np.zeros(len(losses))
for i, v in enumerate(losses):
    pdf[i] += sum(probabilities[values == v])
cdf = np.cumsum(pdf)

i_var = np.argmax(cdf >= 1-alpha)
exact_var = losses[i_var]
exact_cvar = np.dot(pdf[(i_var+1):], losses[(i_var+1):])/sum(pdf[(i_var+1):])

print('Expected Loss E[L]:                %.4f' % expected_loss)
print('Value at Risk VaR[L]:              %.4f' % exact_var)
print('P[L <= VaR[L]]:                    %.4f' % cdf[exact_var])
print('Conditional Value at Risk CVaR[L]: %.4f' % exact_cvar)

Expected Loss E[L]:                0.6409
Value at Risk VaR[L]:              2.0000
P[L <= VaR[L]]:                    0.9591
Conditional Value at Risk CVaR[L]: 3.0000

# plot loss PDF, expected loss, var, and cvar
plt.bar(losses, pdf)
plt.axvline(expected_loss, color='green', linestyle='--', label='E[L]')
plt.axvline(exact_var, color='orange', linestyle='--', label='VaR(L)')
plt.axvline(exact_cvar, color='red', linestyle='--', label='CVaR(L)')
plt.legend(fontsize=15)
plt.xlabel('Loss L ($)', size=15)
plt.ylabel('probability (%)', size=15)
plt.title('Loss Distribution', size=20)
plt.xticks(size=15)
plt.yticks(size=15)
plt.show()

# plot results for Z
plt.plot(z_values, p_z, 'o-', linewidth=3, markersize=8)
plt.grid()
plt.xlabel('Z value', size=15)
plt.ylabel('probability (%)', size=15)
plt.title('Z Distribution', size=20)
plt.xticks(size=15)
plt.yticks(size=15)
plt.show()

# plot results for default probabilities
plt.bar(range(K), p_default)
plt.xlabel('Asset', size=15)
plt.ylabel('probability (%)', size=15)
plt.title('Individual Default Probabilities', size=20)
plt.xticks(range(K), size=15)
plt.yticks(size=15)
plt.grid()
plt.show()

예상 손실

예상 손실을 추정하기 위해 먼저 가중치 합계 연산자를 적용하여 총 손실에 대한 개별 손실을 합산합니다.

\[* mathcal{S}: |x_1, ..., x_K \rangle_K | 0\rangle_{n_S} \mapsto |x_1, ..., x_K \rangle_K | \lambda_1x_1+...+\lambda_K x_K\rangle_{n_S}.\]

결과를 표시하는 데 필요한 필수 비트 수는 다음과 같이 제공됩니다.

\[n_s = \lfloor \log_2( \lambda_1 + ... + \lambda_K ) \rfloor + 1.\]

Once we have the total loss distribution in a quantum register, we can use the techniques described in [Woerner2019] to map a total loss \(L \in \{0, ..., 2^{n_s}-1\}\) to the amplitude of an objective qubit by an operator

\[ | L \rangle_{n_s}|0\rangle \mapsto | L \rangle_{n_s} \left( \sqrt{1 - L/(2^{n_s}-1)}|0\rangle + \sqrt{L/(2^{n_s}-1)}|1\rangle \right),\]

이것은 예상 손실을 평가하기 위해 진폭 추정을 실행할 수 있도록 한다.

# add Z qubits with weight/loss 0
from qiskit.circuit.library import WeightedAdder
agg = WeightedAdder(n_z + K, [0]*n_z + lgd)

from qiskit.circuit.library import LinearAmplitudeFunction

# define linear objective function
breakpoints = [0]
slopes = [1]
offsets = [0]
f_min = 0
f_max = sum(lgd)
c_approx = 0.25

objective = LinearAmplitudeFunction(
    agg.num_sum_qubits,
    slope=slopes,
    offset=offsets,
    # max value that can be reached by the qubit register (will not always be reached)
    domain=(0, 2**agg.num_sum_qubits-1),
    image=(f_min, f_max),
    rescaling_factor=c_approx,
    breakpoints=breakpoints
)

상태 준비 회로를 작성한다.

# define the registers for convenience and readability
qr_state = QuantumRegister(u.num_qubits, 'state')
qr_sum = QuantumRegister(agg.num_sum_qubits, 'sum')
qr_carry = QuantumRegister(agg.num_carry_qubits, 'carry')
qr_obj = QuantumRegister(1, 'objective')

# define the circuit
state_preparation = QuantumCircuit(qr_state, qr_obj, qr_sum, qr_carry, name='A')

# load the random variable
state_preparation.append(u.to_gate(), qr_state)

# aggregate
state_preparation.append(agg.to_gate(), qr_state[:] + qr_sum[:] + qr_carry[:])

# linear objective function
state_preparation.append(objective.to_gate(), qr_sum[:] + qr_obj[:])

# uncompute aggregation
state_preparation.append(agg.to_gate().inverse(), qr_state[:] + qr_sum[:] + qr_carry[:])

# draw the circuit
state_preparation.draw()

             ┌───────┐┌────────┐      ┌───────────┐
    state_0: ┤0      ├┤0       ├──────┤0          ├
             │       ││        │      │           │
    state_1: ┤1      ├┤1       ├──────┤1          ├
             │  P(X) ││        │      │           │
    state_2: ┤2      ├┤2       ├──────┤2          ├
             │       ││        │      │           │
    state_3: ┤3      ├┤3       ├──────┤3          ├
             └───────┘│  adder │┌────┐│  adder_dg │
objective_0: ─────────┤        ├┤2   ├┤           ├
                      │        ││    ││           │
      sum_0: ─────────┤4       ├┤0 F ├┤4          ├
                      │        ││    ││           │
      sum_1: ─────────┤5       ├┤1   ├┤5          ├
                      │        │└────┘│           │
    carry_0: ─────────┤6       ├──────┤6          ├
                      └────────┘      └───────────┘

Before we use QAE to estimate the expected loss, we validate the quantum circuit representing the objective function by just simulating it directly and analyzing the probability of the objective qubit being in the \(|1\rangle\) state, i.e., the value QAE will eventually approximate.

job = execute(state_preparation, backend=Aer.get_backend('statevector_simulator'))

# evaluate resulting statevector
value = 0
for i, a in enumerate(job.result().get_statevector()):
    b = ('{0:0%sb}' % (len(qr_state) + 1)).format(i)[-(len(qr_state) + 1):]
    am = np.round(np.real(a), decimals=4)
    if np.abs(am) > 1e-6 and b[0] == '1':
        value += am**2

print('Exact Expected Loss:   %.4f' % expected_loss)
print('Exact Operator Value:  %.4f' % value)
print('Mapped Operator value: %.4f' % objective.post_processing(value))

Exact Expected Loss:   0.6409
Exact Operator Value:  0.3906
Mapped Operator value: 0.6640

다음은 QAE를 실행하여 고전적인 Monte Carlo 시뮬레이션을 능가하는 2차 속도로 예상 손실을 추정한다.

# set target precision and confidence level
epsilon = 0.01
alpha = 0.05

# construct amplitude estimation
ae = IterativeAmplitudeEstimation(state_preparation=state_preparation,
                                  epsilon=epsilon, alpha=alpha,
                                  objective_qubits=[len(qr_state)],
                                  post_processing=objective.post_processing)
result = ae.run(quantum_instance=Aer.get_backend('qasm_simulator'), shots=100)

# print results
conf_int = np.array(result['confidence_interval'])
print('Exact value:    \t%.4f' % expected_loss)
print('Estimated value:\t%.4f' % result['estimation'])
print('Confidence interval: \t[%.4f, %.4f]' % tuple(conf_int))

Exact value:            0.6409
Estimated value:        0.6913
Confidence interval:    [0.6224, 0.7601]

누적 분포 함수

예상 손실 대신 (고전적 기법을 사용하여 효과적으로 평가할 수 있음) 이제 손실의 누적 분포 함수 (CDF) 를 추정합니다. 이는 기본적으로 기본 자산의 가능한 조합을 모두 평가하거나 Monte Carlo 시뮬레이션에 있는 많은 고전적 샘플을 평가하는 것을 포함한다. QAE에 기반한 알고리즘은 미래에 이러한 분석을 상당히 향상시킬 수 있는 잠재력을 갖고 있다.

To estimate the CDF, i.e., the probability \(\mathbb{P}[L \leq x]\), we again apply \(\mathcal{S}\) to compute the total loss, and then apply a comparator that for a given value \(x\) acts as

\[\begin{split} \mathcal{C}: |L\rangle_n|0> \mapsto \begin{cases} |L\rangle_n|1> & \text{if}\quad L \leq x \\ |L\rangle_n|0> & \text{if}\quad L > x. \end{cases}\end{split}\]

결과적인 양자 상태는 다음과 같이 기록될 수 있다.

\[ \sum_{L = 0}^{x} \sqrt{p_{L}}|L\rangle_{n_s}|1\rangle + \sum_{L = x+1}^{2^{n_s}-1} \sqrt{p_{L}}|L\rangle_{n_s}|1\rangle,\]

여기서 우리는 불확실성 모델의 세부사항들을 제시하는 대신에 합산된 손실 값들 및 대응하는 확률들을 추정한다.

The CDF(\(x\)) equals the probability of measuring \(|1\rangle\) in the objective qubit and QAE can be directly used to estimate it.

# set x value to estimate the CDF
x_eval = 2

comparator = IntegerComparator(agg.num_sum_qubits, x_eval + 1, geq=False)
comparator.draw()

  state_0: ──■──────────────■──
             │              │
  state_1: ──┼────■─────────┼──
             │  ┌─┴─┐┌───┐  │
compare_0: ──┼──┤ X ├┤ X ├──┼──
           ┌─┴─┐└─┬─┘└───┘┌─┴─┐
     a0_0: ┤ X ├──■───────┤ X ├
           └───┘          └───┘

def get_cdf_circuit(x_eval):
    # define the registers for convenience and readability
    qr_state = QuantumRegister(u.num_qubits, 'state')
    qr_sum = QuantumRegister(agg.num_sum_qubits, 'sum')
    qr_carry = QuantumRegister(agg.num_carry_qubits, 'carry')
    qr_obj = QuantumRegister(1, 'objective')
    qr_compare = QuantumRegister(1, 'compare')

    # define the circuit
    state_preparation = QuantumCircuit(qr_state, qr_obj, qr_sum, qr_carry, name='A')

    # load the random variable
    state_preparation.append(u, qr_state)

    # aggregate
    state_preparation.append(agg, qr_state[:] + qr_sum[:] + qr_carry[:])

    # comparator objective function
    comparator = IntegerComparator(agg.num_sum_qubits, x_eval + 1, geq=False)
    state_preparation.append(comparator, qr_sum[:] + qr_obj[:] + qr_carry[:])

    # uncompute aggregation
    state_preparation.append(agg.inverse(), qr_state[:] + qr_sum[:] + qr_carry[:])

    return state_preparation

state_preparation = get_cdf_circuit(x_eval)

다시, 양자회로의 유효성을 검증하기 위해 양자 시뮬레이션을 먼저 사용합니다.

job = execute(state_preparation, backend=Aer.get_backend('statevector_simulator'))

state_preparation.draw()

             ┌───────┐┌────────┐        ┌───────────┐
    state_0: ┤0      ├┤0       ├────────┤0          ├
             │       ││        │        │           │
    state_1: ┤1      ├┤1       ├────────┤1          ├
             │  P(X) ││        │        │           │
    state_2: ┤2      ├┤2       ├────────┤2          ├
             │       ││        │        │           │
    state_3: ┤3      ├┤3       ├────────┤3          ├
             └───────┘│  adder │┌──────┐│  adder_dg │
objective_0: ─────────┤        ├┤2     ├┤           ├
                      │        ││      ││           │
      sum_0: ─────────┤4       ├┤0     ├┤4          ├
                      │        ││  cmp ││           │
      sum_1: ─────────┤5       ├┤1     ├┤5          ├
                      │        ││      ││           │
    carry_0: ─────────┤6       ├┤3     ├┤6          ├
                      └────────┘└──────┘└───────────┘

# evaluate resulting statevector
var_prob = 0
for i, a in enumerate(job.result().get_statevector()):
    b = ('{0:0%sb}' % (len(qr_state) + 1)).format(i)[-(len(qr_state) + 1):]
    prob = np.abs(a)**2
    if prob > 1e-6 and b[0] == '1':
        var_prob += prob
print('Operator CDF(%s)' % x_eval + ' = %.4f' % var_prob)
print('Exact    CDF(%s)' % x_eval + ' = %.4f' % cdf[x_eval])

Operator CDF(2) = 0.9591
Exact    CDF(2) = 0.9591

그 다음에 QAE를 실행하여 주어진 :math:’x’ 에 대한 CDF를 추정한다.

# set target precision and confidence level
epsilon = 0.01
alpha = 0.05

# construct amplitude estimation
ae_cdf = IterativeAmplitudeEstimation(state_preparation=state_preparation,
                                      epsilon=epsilon, alpha=alpha,
                                      objective_qubits=[len(qr_state)])
result_cdf = ae_cdf.run(quantum_instance=Aer.get_backend('qasm_simulator'), shots=100)

# print results
conf_int = np.array(result_cdf['confidence_interval'])
print('Exact value:    \t%.4f' % cdf[x_eval])
print('Estimated value:\t%.4f' % result_cdf['estimation'])
print('Confidence interval: \t[%.4f, %.4f]' % tuple(conf_int))

Exact value:            0.9591
Estimated value:        0.9595
Confidence interval:    [0.9584, 0.9605]

최대예상손실액 (Value at Risk)

다음에는 이중 검색 및 QAE를 사용하여 위험에 있는 값을 추정하기 위해 CDF를 효율적으로 평가합니다.

def run_ae_for_cdf(x_eval, epsilon=0.01, alpha=0.05, simulator='qasm_simulator'):

    # construct amplitude estimation
    state_preparation = get_cdf_circuit(x_eval)
    ae_var = IterativeAmplitudeEstimation(state_preparation=state_preparation,
                                          epsilon=epsilon, alpha=alpha,
                                          objective_qubits=[len(qr_state)])
    result_var = ae_var.run(quantum_instance=Aer.get_backend(simulator), shots=100)

    return result_var['estimation']

def bisection_search(objective, target_value, low_level, high_level, low_value=None, high_value=None):
    """
    Determines the smallest level such that the objective value is still larger than the target
    :param objective: objective function
    :param target: target value
    :param low_level: lowest level to be considered
    :param high_level: highest level to be considered
    :param low_value: value of lowest level (will be evaluated if set to None)
    :param high_value: value of highest level (will be evaluated if set to None)
    :return: dictionary with level, value, num_eval
    """

    # check whether low and high values are given and evaluated them otherwise
    print('--------------------------------------------------------------------')
    print('start bisection search for target value %.3f' % target_value)
    print('--------------------------------------------------------------------')
    num_eval = 0
    if low_value is None:
        low_value = objective(low_level)
        num_eval += 1
    if high_value is None:
        high_value = objective(high_level)
        num_eval += 1

    # check if low_value already satisfies the condition
    if low_value > target_value:
        return {'level': low_level, 'value': low_value, 'num_eval': num_eval, 'comment': 'returned low value'}
    elif low_value == target_value:
        return {'level': low_level, 'value': low_value, 'num_eval': num_eval, 'comment': 'success'}

    # check if high_value is above target
    if high_value < target_value:
        return {'level': high_level, 'value': high_value, 'num_eval': num_eval, 'comment': 'returned low value'}
    elif high_value == target_value:
        return {'level': high_level, 'value': high_value, 'num_eval': num_eval, 'comment': 'success'}

    # perform bisection search until
    print('low_level    low_value    level    value    high_level    high_value')
    print('--------------------------------------------------------------------')
    while high_level - low_level > 1:

        level = int(np.round((high_level + low_level) / 2.0))
        num_eval += 1
        value = objective(level)

        print('%2d           %.3f        %2d       %.3f    %2d            %.3f' \
              % (low_level, low_value, level, value, high_level, high_value))

        if value >= target_value:
            high_level = level
            high_value = value
        else:
            low_level = level
            low_value = value

    # return high value after bisection search
    print('--------------------------------------------------------------------')
    print('finished bisection search')
    print('--------------------------------------------------------------------')
    return {'level': high_level, 'value': high_value, 'num_eval': num_eval, 'comment': 'success'}

# run bisection search to determine VaR
objective = lambda x: run_ae_for_cdf(x)
bisection_result = bisection_search(objective, 1-alpha, min(losses)-1, max(losses), low_value=0, high_value=1)
var = bisection_result['level']

--------------------------------------------------------------------
start bisection search for target value 0.950
--------------------------------------------------------------------
low_level    low_value    level    value    high_level    high_value
--------------------------------------------------------------------
-1           0.000         1       0.751     3            1.000
 1           0.751         2       0.960     3            1.000
--------------------------------------------------------------------
finished bisection search
--------------------------------------------------------------------

print('Estimated Value at Risk: %2d' % var)
print('Exact Value at Risk:     %2d' % exact_var)
print('Estimated Probability:    %.3f' % bisection_result['value'])
print('Exact Probability:        %.3f' % cdf[exact_var])

Estimated Value at Risk:  2
Exact Value at Risk:      2
Estimated Probability:    0.960
Exact Probability:        0.959

위험에 있는 조건부 값

Last, we compute the CVaR, i.e. the expected value of the loss conditional to it being larger than or equal to the VaR. To do so, we evaluate a piecewise linear objective function \(f(L)\), dependent on the total loss \(L\), that is given by

\[\begin{split} f(L) = \begin{cases} 0 & \text{if}\quad L \leq VaR \\ L & \text{if}\quad L > VaR. \end{cases}\end{split}\]

To normalize, we have to divide the resulting expected value by the VaR-probability, i.e. \(\mathbb{P}[L \leq VaR]\).

# define linear objective
breakpoints = [0, var]
slopes = [0, 1]
offsets = [0, 0]  # subtract VaR and add it later to the estimate
f_min = 0
f_max = 3 - var
c_approx = 0.25

cvar_objective = LinearAmplitudeFunction(
    agg.num_sum_qubits,
    slopes,
    offsets,
    domain=(0, 2**agg.num_sum_qubits - 1),
    image=(f_min, f_max),
    rescaling_factor=c_approx,
    breakpoints=breakpoints
)

cvar_objective.draw()

        ┌─────────┐┌──────┐┌─────────┐┌─────────┐
q138_0: ┤0        ├┤0     ├┤0        ├┤0        ├
        │         ││      ││         ││         │
q138_1: ┤1 LinRot ├┤1     ├┤1 LinRot ├┤1        ├
        │         ││      ││         ││         │
q139_0: ┤2        ├┤  cmp ├┤2        ├┤  cmp_dg ├
        └─────────┘│      │└────┬────┘│         │
  a4_0: ───────────┤2     ├─────■─────┤2        ├
                   │      │           │         │
  a4_1: ───────────┤3     ├───────────┤3        ├
                   └──────┘           └─────────┘

# define the registers for convenience and readability
qr_state = QuantumRegister(u.num_qubits, 'state')
qr_sum = QuantumRegister(agg.num_sum_qubits, 'sum')
qr_carry = QuantumRegister(agg.num_carry_qubits, 'carry')
qr_obj = QuantumRegister(1, 'objective')
qr_work = QuantumRegister(cvar_objective.num_ancillas - len(qr_carry), 'work')

# define the circuit
state_preparation = QuantumCircuit(qr_state, qr_obj, qr_sum, qr_carry, qr_work, name='A')

# load the random variable
state_preparation.append(u, qr_state)

# aggregate
state_preparation.append(agg, qr_state[:] + qr_sum[:] + qr_carry[:])

# linear objective function
state_preparation.append(cvar_objective, qr_sum[:] + qr_obj[:] + qr_carry[:] + qr_work[:])

# uncompute aggregation
state_preparation.append(agg.inverse(), qr_state[:] + qr_sum[:] + qr_carry[:])

<qiskit.circuit.instructionset.InstructionSet at 0x7fdc6a68f110>

다시, 양자회로의 유효성을 검증하기 위해 양자 시뮬레이션을 먼저 사용합니다.

job = execute(state_preparation, backend=Aer.get_backend('statevector_simulator'))

# evaluate resulting statevector
value = 0
for i, a in enumerate(job.result().get_statevector()):
    b = ('{0:0%sb}' % (len(qr_state) + 1)).format(i)[-(len(qr_state) + 1):]
    am = np.round(np.real(a), decimals=4)
    if np.abs(am) > 1e-6 and b[0] == '1':
        value += am**2

# normalize and add VaR to estimate
value = cvar_objective.post_processing(value)
d = (1.0 - bisection_result['value'])
v = value / d if d != 0 else 0
normalized_value = v + var
print('Estimated CVaR: %.4f' % normalized_value)
print('Exact CVaR:     %.4f' % exact_cvar)

Estimated CVaR: 3.3121
Exact CVaR:     3.0000

다음에는 QAE를 실행하여 CVS를 추정한다.

# set target precision and confidence level
epsilon = 0.01
alpha = 0.05

# construct amplitude estimation
ae_cvar = IterativeAmplitudeEstimation(state_preparation=state_preparation,
                                       epsilon=epsilon, alpha=alpha,
                                       objective_qubits=[len(qr_state)],
                                       post_processing=cvar_objective.post_processing)
result_cvar = ae_cvar.run(quantum_instance=Aer.get_backend('qasm_simulator'), shots=100)

# print results
d = (1.0 - bisection_result['value'])
v = result_cvar['estimation'] / d if d != 0 else 0
print('Exact CVaR:    \t%.4f' % exact_cvar)
print('Estimated CVaR:\t%.4f' % (v + var))

Exact CVaR:     3.0000
Estimated CVaR: 3.2474

import qiskit.tools.jupyter
%qiskit_version_table
%qiskit_copyright

Version Information

Qiskit Software	Version
Qiskit	None
Terra	0.17.0.dev0+4ada179
Aer	0.6.1
Ignis	0.5.0.dev0+470d8cc
Aqua	0.9.0.dev0+5a88b59
IBM Q Provider	0.8.0
System information
Python	3.7.7 (default, May 6 2020, 04:59:01) [Clang 4.0.1 (tags/RELEASE_401/final)]
OS	Darwin
CPUs	2
Memory (Gb)	16.0
Tue Oct 20 11:03:45 2020 CEST

This code is a part of Qiskit

© Copyright IBM 2017, 2020.

This code is licensed under the Apache License, Version 2.0. You may
obtain a copy of this license in the LICENSE.txt file in the root directory
of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.

Any modifications or derivative works of this code must retain this
copyright notice, and modified files need to carry a notice indicating
that they have been altered from the originals.