[from qiskit] Error in MaximumLikelihoodAmplitudeEstimation when using Sampler

[ElePT]

What should we add?

Environment

• Qiskit Terra version: 0.22.2
• Python version: 3.9.13
• Operating system: macOS

What is happening?

Wrong answer when using MaximumLikelihoodAmplitudeEstimationwith with Sampler

How can we reproduce the issue?

from qiskit.utils import QuantumInstance
from qiskit.algorithms import EstimationProblem, MaximumLikelihoodAmplitudeEstimation
from qiskit import Aer, QuantumCircuit
from qiskit.primitives import Sampler
import numpy as np

number_qubits = 3
state_prep =  QuantumCircuit(number_qubits)
for i in range(number_qubits):
    state_prep.h(i)

qi = QuantumInstance(Aer.get_backend("aer_simulator"), shots=4000)
sampler = Sampler()

problem = EstimationProblem(
    state_preparation=state_prep,
    objective_qubits=[number_qubits-1],
)

ae_instance = MaximumLikelihoodAmplitudeEstimation(evaluation_schedule=[0,1,2], quantum_instance=qi)
ae_sampler = MaximumLikelihoodAmplitudeEstimation(evaluation_schedule=[0,1,2], sampler=sampler)

result_instance  = ae_instance.estimate(problem)
result_sampler = ae_sampler.estimate(problem)

print("--------- Correct result ---------")
print("Good counts with aer simulator: ",result_instance.good_counts)
print("Circuit results with aer simulator: " ,result_instance.circuit_results)

conf_int = np.array(result_instance.confidence_interval)
print("Estimated value with aer simulator:    \t%.4f" % (result_instance.estimation_processed))
print("Confidence interval with aer simulator :\t[%.4f, %.4f]" % tuple(conf_int))

print("\n")

print("--------- Incorrect result ---------")
print("Good counts with Sampler:",result_sampler.good_counts)
print("Circuit results with Sampler:",result_sampler.circuit_results)

conf_int = np.array(result_sampler.confidence_interval)
print("Estimated value with Sampler:    \t%.4f" % (result_sampler.estimation_processed))
print("Confidence interval with Sampler:\t[%.4f, %.4f]" % tuple(conf_int))

--------- Correct result ---------
Good counts with aer simulator:  [1976, 2000, 2061]
Circuit results with aer simulator:  
[{'0': 2024, '1': 1976}, 
{'0': 2000, '1': 2000},
 {'0': 1939, '1': 2061}]
Estimated value with aer simulator:  0.5020
Confidence interval with aer simulator : [0.4994, 0.5046]

--------- Incorrect result ---------
Good counts with Sampler: [0, 0, 0]
Circuit results with Sampler: 
[{'000': 0.4999999999999998, '001': 0.4999999999999998}, 
{'000': 0.49999999999999895, '001': 0.49999999999999895},
 {'000': 0.49999999999999817, '001': 0.49999999999999817}]
Estimated value with Sampler:  0.0000
Confidence interval with Sampler:[-0.0000, 0.0000]

What should happen?

Using Sampler should get similar results as with aer simulator

Any suggestions?

I have identified the problem in MaximumLikelihoodAmplitudeEstimation and have created a diagram explaining it:

To fix the issue, I suggest changing some parts in the function estimate in MaximumLikelihoodAmplitudeEstimation. From this:

                if shots is None:
                    for i, quasi_dist in enumerate(ret.quasi_dists):
                        circuit_result = {
                            np.binary_repr(k, circuits[i].num_qubits): v
                            for k, v in quasi_dist.items()
                        }
                        result.circuit_results.append(circuit_result)
                    shots = 1
                else:
                    # get counts and construct MLE input
                    for circuit in circuits:
                        counts = {
                            np.binary_repr(k, circuit.num_qubits): round(v * shots)
                            for k, v in ret.quasi_dists[0].items()
                        }
                        result.circuit_results.append(counts)

To this:

                if shots is None:
                    for i, quasi_dist in enumerate(ret.quasi_dists):
                        circuit_result = {
                            str(k): v
                            for k, v in quasi_dist.items()
                        }
                        result.circuit_results.append(circuit_result)
                    shots = 1
                else:
                    # get counts and construct MLE input
                    for index,circuit in enumerate(circuits):
                        counts = {
                            str(k): round(v * shots)
                            for k, v in ret.quasi_dists[index].items()
                        }
                        result.circuit_results.append(counts)

If everyone agrees with the solution, I would love to try and push the changes (it will be my first time contributing to Qiskit)