From 485133cd843d5bb46dcb112c8ede38d520b5dcc2 Mon Sep 17 00:00:00 2001 From: Manan Gupta Date: Sun, 9 Aug 2026 12:00:46 -0700 Subject: [PATCH] Fix QPCA crash on 4-feature input and k=2 returning 0 --- pyqpanda-algorithm/pyqpanda_alg/QPCA/QPCA.py | 104 +++++++++---------- test/QPCA/Test_qpca.py | 73 +++++++++++++ 2 files changed, 125 insertions(+), 52 deletions(-) diff --git a/pyqpanda-algorithm/pyqpanda_alg/QPCA/QPCA.py b/pyqpanda-algorithm/pyqpanda_alg/QPCA/QPCA.py index 1a485862..1f45aaeb 100644 --- a/pyqpanda-algorithm/pyqpanda_alg/QPCA/QPCA.py +++ b/pyqpanda-algorithm/pyqpanda_alg/QPCA/QPCA.py @@ -13,7 +13,6 @@ from pyqpanda3.core import CPUQVM, QCircuit, QProg, TOFFOLI, SWAP, CNOT, U1, U3, CR, I, H, X, Y, measure import numpy as np import math -import sys from .. plugin import * @@ -201,15 +200,22 @@ def qpca(sample_A, k): """ QPCA is a quantum version of the classical PCA algorithm, which is widely used in data analysis and machine learning. + Only 2-feature input is supported. The n=4 branches of the circuit helpers + are incomplete, so a larger input raises NotImplementedError. + Parameters: sample_A: ``ndarray``\n - the input matrix for analysis + the input matrix for analysis, shape (n_samples, 2) k: ``int``\n - the dimension to reduce + the dimension to reduce, 1 or 2 Returns: out: ``ndarray``\n - the output matrix after reducing dimension + the output matrix after reducing dimension, shape (n_samples, k) + + Raises: + ValueError: k is not 1 or 2.\n + NotImplementedError: sample_A has other than 2 features. Examples: .. code-block:: python @@ -222,61 +228,55 @@ def qpca(sample_A, k): print(data_q) """ norm_x, A = _preprocessing(sample_A) - lambda_A, vector_A = np.linalg.eig(A) - A1 = A.reshape(1, 2 ** A.shape[0])[0] + n = A.shape[0] + if n != 2: + raise NotImplementedError( + 'qpca supports 2-feature input only, got %d features' % n) + if k not in (1, 2): + raise ValueError('k must be 1 or 2, got %r' % (k,)) + + # the covariance matrix is symmetric, so eigh gives real eigenvalues and an + # orthonormal basis. sorted descending, principal component first. + lambda_A, vector_A = np.linalg.eigh(A) + order = np.argsort(lambda_A)[::-1] + lambda_A = lambda_A[order] + vector_A = vector_A[:, order] + A1 = A.reshape(1, 2 ** n)[0] state_vector = np.sqrt(1 / np.sum(A1 * A1)) * A1 - if A.shape[0] == 2: - if k == 2: - tao = min(lambda_A) - 1 - elif k == 1: - tao = np.sum(lambda_A)/2 - else: - print('The K Error!') - sys.exit() - qm = _QMachine(5, 5) - if A.shape[0] == 4: - qm = _QMachine(8, 8) + # tao is the eigenvalue threshold. k=1 keeps the top eigenvalue only, + # k=2 sits below both and keeps the whole space. + if k == 2: + tao = min(lambda_A) - 1 + else: + tao = np.sum(lambda_A) / 2 + qm = _QMachine(5, 5) prog = QProg() cir = QCircuit() - _init_cir(qm, state_vector, A.shape[0]) - cir << _phase_estimation_cir(qm.q_list, lambda_A, tao, A.shape[0]) - cir << _transition_cir(qm.q_list, A.shape[0]) + _init_cir(qm, state_vector, n) + cir << _phase_estimation_cir(qm.q_list, lambda_A, tao, n) + cir << _transition_cir(qm.q_list, n) cir << _cnot_cir(qm.q_list) - cir << _transition_reverse_cir(qm.q_list, A.shape[0]) - cir << _phase_estimation_reverse_cir(qm.q_list, lambda_A, tao, A.shape[0]) + cir << _transition_reverse_cir(qm.q_list, n) + cir << _phase_estimation_reverse_cir(qm.q_list, lambda_A, tao, n) prog << cir prog <<_measure_cir(prog, qm.q_list, qm.q_list) - # result = qm.machine.run_with_configuration(prog, qm.c_list, 8192) qm.machine.run(prog, 8192) result = qm.machine.result().get_prob_dict(qm.q_list) - a = [] - data = 0 - if A.shape[0] == 2: - for i, v in enumerate(result.keys()): - if int(v[-1]) == 1: - # a.append(float("%.4f" % (result[v] / 8192))) - a.append(float("%.4f" % (result[v]))) - if k == 2: - result_idealA = state_vector - if k == 1: - i = np.argmax(lambda_A) - result_idealA = (np.kron(vector_A[:, i], vector_A[:, i]) * lambda_A[i]) / ( - np.sqrt(lambda_A[i] * lambda_A[i])) - result_idealA1 = [] - for i in range(len(result_idealA)): - if result_idealA[i] != 0: - result_idealA1.append(result_idealA[i]) - result_idealA = result_idealA1 - sum_A = np.sum(np.array(a)) - result_circuitA = [] - for i in range(len(a)): - result_circuitA.append(np.sqrt(a[i] / sum_A)) - if k == 1: - vector_qpca = [] - vector_qpca.append(np.sqrt(result_circuitA[0])) - vector_qpca.append(np.sqrt(result_circuitA[-1])) - vector_qpca = np.array(vector_qpca).reshape(1, 2) - data = np.dot(norm_x, np.transpose(vector_qpca)) - return data + + if k == 2: + # tao is below every eigenvalue, so post-selection succeeds with + # probability 1 and nothing is filtered out. the retained subspace is + # the whole feature space and its basis comes from vector_A. + return np.dot(norm_x, vector_A) + + # post-select the ancilla, qubit 0, which is the last character of the key. + # sorted() so the branch order does not depend on dict insertion order. + a = [result[v] for v in sorted(result) if int(v[-1]) == 1] + sum_A = np.sum(np.array(a)) + # probabilities back to amplitudes. unit norm when post-selection leaves + # two outcomes, which is the non-degenerate case. + result_circuitA = [np.sqrt(x / sum_A) for x in a] + vector_qpca = np.array([result_circuitA[0], result_circuitA[-1]]).reshape(1, 2) + return np.dot(norm_x, np.transpose(vector_qpca)) diff --git a/test/QPCA/Test_qpca.py b/test/QPCA/Test_qpca.py index ca2c40d3..0a757d8f 100644 --- a/test/QPCA/Test_qpca.py +++ b/test/QPCA/Test_qpca.py @@ -8,6 +8,15 @@ sys.path.append((Path.cwd().parent.parent).__str__()) +def _classical_pca(x): + """Reference PCA: eigendecomposition of the covariance matrix, descending.""" + norm_x = x - x.mean(axis=0) + cov = np.dot(norm_x.T, norm_x) + w, v = np.linalg.eigh(cov) + order = np.argsort(w)[::-1] + return norm_x, w[order], v[:, order] + + class TestQPCA: """QPCA测试类""" @@ -16,6 +25,11 @@ def sample_data(self): """提供标准测试数据""" return np.array([[-1, 2], [-2, -1], [-1, -2], [1, 3], [2, 1], [3, 2]]) + @pytest.fixture + def sample_data_4d(self): + rng = np.random.RandomState(0) + return rng.randn(20, 4) + def test_qpca_normal(self, sample_data): """测试QPCA正常功能""" from pyqpanda_alg.QPCA import qpca @@ -26,5 +40,64 @@ def test_qpca_normal(self, sample_data): assert data_q.shape == (6,1) assert isinstance(data_q, np.ndarray) + def test_qpca_k1_matches_classical_pca(self, sample_data): + from pyqpanda_alg.QPCA import qpca + + data_q = qpca(sample_data, 1)[:, 0] + norm_x, _, v = _classical_pca(sample_data) + expected = np.dot(norm_x, v[:, 0]) + if np.dot(data_q, expected) < 0: + expected = -expected + + # sampling noise on 8192 shots keeps this loose + assert np.corrcoef(data_q, expected)[0, 1] > 0.999 + assert np.max(np.abs(data_q - expected)) < 0.15 + + def test_qpca_k1_direction_is_unit_norm(self, sample_data): + """Holds for this fixture, not in general: a degenerate spectrum leaves + the post-selection with a single outcome and the norm comes out at + sqrt(2). Pinned to the documented example on purpose.""" + from pyqpanda_alg.QPCA import qpca + + norm_x, _, _ = _classical_pca(sample_data) + data_q = qpca(sample_data, 1) + # recover the projection direction the circuit produced + vec = np.linalg.lstsq(norm_x, data_q, rcond=None)[0][:, 0] + assert np.linalg.norm(vec) == pytest.approx(1.0, abs=1e-6) + + def test_qpca_k2_returns_pca_scores(self, sample_data): + from pyqpanda_alg.QPCA import qpca + + data_q = qpca(sample_data, 2) + norm_x, w, v = _classical_pca(sample_data) + + assert isinstance(data_q, np.ndarray) + assert data_q.shape == (6, 2) + + # PCA scores are decorrelated and carry the eigenvalues as their + # column variances, in descending order. neither holds for the + # centred input, so this pins the basis and not just the subspace. + cov_q = np.dot(data_q.T, data_q) + assert abs(cov_q[0, 1]) < 1e-9 + assert np.allclose(np.diag(cov_q), w) + + expected = np.dot(norm_x, v) + for j in range(2): + sign = np.sign(np.dot(data_q[:, j], expected[:, j])) + assert np.allclose(data_q[:, j], sign * expected[:, j]) + + @pytest.mark.parametrize('k', [0, 3, -1]) + def test_qpca_invalid_k_raises(self, sample_data, k): + from pyqpanda_alg.QPCA import qpca + + with pytest.raises(ValueError): + qpca(sample_data, k) + + def test_qpca_four_features_raises(self, sample_data_4d): + from pyqpanda_alg.QPCA import qpca + with pytest.raises(NotImplementedError): + qpca(sample_data_4d, 1) + with pytest.raises(NotImplementedError): + qpca(sample_data_4d, 2)