Qiskit Runtime Primitives: Sampler and Estimator

Why Primitives?

For several years, Qiskit code that ran circuits on IBM backends used execute() or the now-deprecated backend.run() pattern. That API was informal: you passed a transpiled circuit, called .result(), and parsed a dictionary of counts. It worked, but it created a problem. The same code could not move cleanly between a local simulator and real hardware without manually wrangling result formats, shots, and transpilation steps.

Qiskit Runtime primitives replace that pattern with two well-defined interfaces:

• Sampler gets bitstring distributions. You give it circuits and it returns quasi-probability distributions over measurement outcomes.
• Estimator gets expectation values. You give it circuits, observables, and parameter values, and it returns <ψ|O|ψ> with standard error estimates.

Both primitives run identically against local statevector simulators and against real IBM Quantum hardware. The only change between environments is which class you import.

Architecture Deep Dive: What “Primitive” Means

The word “primitive” here carries a specific architectural meaning. A primitive is the smallest, most fundamental operation a quantum computing service exposes. Instead of offering a low-level “run this circuit on that backend” endpoint, the runtime exposes exactly two primitives: sample and estimate. Every quantum workload decomposes into one of those two tasks.

The primitives provide a hardware-agnostic abstraction layer. Your code describes what you want to compute (bitstring samples or expectation values), and the primitive handles how to execute that computation on a specific backend. This separation means the same user code runs on a statevector simulator, a noisy simulator, or a 127-qubit Eagle processor without modification beyond the import line.

Execution Model

Under the old backend.run() API, the caller was responsible for transpiling circuits, managing shots, submitting jobs, and parsing result dictionaries. The primitives shift that responsibility. The execution flow looks like this:

 ┌────────────┐
 │ User Code  │   You build circuits, observables, parameters
 └─────┬──────┘
       │  list of PUBs
       ▼
 ┌────────────┐
 │ Primitive   │   Sampler or Estimator validates inputs,
 │ (V2 API)   │   handles transpilation and qubit mapping
 └─────┬──────┘
       │  ISA circuits (native gates, physical qubits)
       ▼
 ┌────────────┐
 │  Backend   │   QPU or simulator executes the circuits
 └─────┬──────┘
       │  raw measurement data
       ▼
 ┌────────────┐
 │  Results   │   BitArray (Sampler) or evs/stds (Estimator)
 └────────────┘

Key differences from the old API:

1. Transpilation is handled internally when using the runtime primitives on real hardware. You can still pre-transpile for control, but the primitive validates and re-transpiles if needed.
2. Qubit mapping is automatic. The primitive knows which physical qubits your logical qubits map to.
3. Job submission is batched. Multiple PUBs go in one .run() call, producing one job on the queue rather than many.
4. Result extraction is typed. You get a BitArray or structured expectation values, not a raw dictionary that differs between simulator and hardware.

Installing and Importing

pip install qiskit qiskit-ibm-runtime

For local simulation, Qiskit provides statevector-backed primitives with no account required:

from qiskit.primitives import StatevectorSampler, StatevectorEstimator

For real hardware, use the runtime versions:

from qiskit_ibm_runtime import QiskitRuntimeService, SamplerV2, EstimatorV2

Saving Your IBM Quantum Account

To access real hardware, save your API token once:

from qiskit_ibm_runtime import QiskitRuntimeService

QiskitRuntimeService.save_account(
    channel="ibm_quantum",
    token="YOUR_IBM_QUANTUM_TOKEN",   # from quantum.ibm.com
    overwrite=True,
)

On subsequent runs, load with:

service = QiskitRuntimeService()
backend = service.least_busy(operational=True, simulator=False)

The PUB Format: Complete Reference

Both primitives accept inputs as PUBs: Primitive Unified Blocs. A PUB is a tuple containing the inputs for a single execution unit. You pass a list of PUBs to the primitive so you can batch multiple jobs in one call.

Sampler PUB Formats

The Sampler accepts three PUB formats:

from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.primitives import StatevectorSampler

sampler = StatevectorSampler()

# --- Format 1: (circuit,) ---
# No parameters, default shots
qc = QuantumCircuit(2)
qc.h(0)
qc.cx(0, 1)
qc.measure_all()

job = sampler.run([(qc,)])
result = job.result()
print(result[0].data.meas.get_counts())

# --- Format 2: (circuit, parameter_values) ---
# Parameterized circuit with bound values
theta = ParameterVector("t", 2)
pqc = QuantumCircuit(2)
pqc.ry(theta[0], 0)
pqc.ry(theta[1], 1)
pqc.measure_all()

# parameter_values shape: (num_sets, num_params)
# Single parameter set: [[0.5, 1.2]]
job = sampler.run([(pqc, [[0.5, 1.2]])])
result = job.result()
print(result[0].data.meas.get_counts())

# --- Format 3: (circuit, parameter_values, shots) ---
# Explicit shot count
job = sampler.run([(pqc, [[0.5, 1.2]], 4096)])
result = job.result()
print(result[0].data.meas.get_counts())

Estimator PUB Formats

The Estimator accepts four PUB formats:

import numpy as np
from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import StatevectorEstimator

estimator = StatevectorEstimator()

theta = ParameterVector("t", 1)
qc = QuantumCircuit(1)
qc.ry(theta[0], 0)

obs = SparsePauliOp("Z")

# --- Format 1: (circuit, observables) ---
# Circuit must have no free parameters
qc_bound = QuantumCircuit(1)
qc_bound.ry(np.pi / 4, 0)
job = estimator.run([(qc_bound, obs)])
result = job.result()
print(result[0].data.evs)  # scalar expectation value

# --- Format 2: (circuit, observables, parameter_values) ---
job = estimator.run([(qc, obs, [[np.pi / 4]])])
result = job.result()
print(result[0].data.evs)

# --- Format 3: (circuit, observables, parameter_values, precision) ---
# precision controls the target standard error
job = estimator.run([(qc, obs, [[np.pi / 4]], 0.01)])
result = job.result()
print(result[0].data.evs, result[0].data.stds)

Parameter Values Shape

The parameter_values array shape must be (num_parameter_sets, num_parameters). If your circuit has a ParameterVector of length 3, each inner list must have exactly 3 elements:

# ParameterVector of length 3
params = ParameterVector("p", 3)
qc = QuantumCircuit(2)
qc.ry(params[0], 0)
qc.ry(params[1], 1)
qc.rz(params[2], 0)
qc.measure_all()

# Correct: each inner list has 3 values matching 3 parameters
sampler.run([(qc, [[0.1, 0.2, 0.3]])])            # 1 parameter set
sampler.run([(qc, [[0.1, 0.2, 0.3],
                    [0.4, 0.5, 0.6]])])             # 2 parameter sets

Sampler: Bitstring Distributions

The Sampler executes circuits that end in measurements and returns quasi-probability distributions over the resulting bitstrings.

Bell State Example

from qiskit import QuantumCircuit
from qiskit.primitives import StatevectorSampler

# Build a Bell state circuit
qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

# Create the Sampler and run
sampler = StatevectorSampler()
job = sampler.run([(qc, [], 1024)])   # one PUB: circuit, no params, 1024 shots
result = job.result()

# Access the BitArray for the first PUB
pub_result = result[0]
bit_array = pub_result.data.c    # 'c' matches the ClassicalRegister name

# Convert to a counts dictionary
counts = bit_array.get_counts()
print(counts)
# Expected: {'00': ~512, '11': ~512}

# Or get quasi-probabilities directly
quasi_probs = pub_result.data.c.get_int_counts()
print(quasi_probs)
# {0: ~512, 3: ~512}  (binary 00 and 11)

The result is a BitArray, not a raw dictionary. This is intentional: the BitArray tracks the full shot-level data and lets you compute counts, probabilities, or marginals over any subset of bits.

BitArray Deep Dive

The BitArray is the core result type for the Sampler. It stores every individual shot as a bitstring and provides methods to analyze the data in different ways.

BitArray API Reference

from qiskit import QuantumCircuit
from qiskit.primitives import StatevectorSampler

qc = QuantumCircuit(3)
qc.h(0)
qc.cx(0, 1)
qc.cx(0, 2)
qc.measure_all()

sampler = StatevectorSampler()
result = sampler.run([(qc, [], 2048)]).result()
ba = result[0].data.meas

# .get_counts() returns {bitstring: count}
print(ba.get_counts())
# {'000': ~1024, '111': ~1024}

# .get_int_counts() returns {integer: count}
print(ba.get_int_counts())
# {0: ~1024, 7: ~1024}

# .get_bitstrings() returns a list of bitstrings, one per shot
bitstrings = ba.get_bitstrings()
print(len(bitstrings))   # 2048
print(bitstrings[0])     # e.g. '000' or '111'

# .num_bits gives the number of classical bits
print(ba.num_bits)  # 3

# .num_shots gives the number of shots
print(ba.num_shots)  # 2048

Marginal Distributions with Slicing

You can compute marginal distributions over a subset of qubits. For example, to look at only qubit 0 and ignore qubits 1 and 2:

# Slice the BitArray to keep only specific bit indices
# bit index 0 is the least significant bit
marginal_ba = ba.slice_bits([0])
print(marginal_ba.get_counts())
# {'0': ~1024, '1': ~1024}

# Keep qubits 0 and 2 only (skip qubit 1)
marginal_02 = ba.slice_bits([0, 2])
print(marginal_02.get_counts())
# {'00': ~1024, '11': ~1024}

Shot-Level Boolean Arrays for Post-Processing

For classical post-processing, you can convert the BitArray into a NumPy array of integers and perform bitwise operations:

import numpy as np

# Get all bitstrings as a list
bitstrings = ba.get_bitstrings()

# Convert to a boolean array: did qubit 0 measure |1>?
qubit0_results = np.array([int(bs[-1]) for bs in bitstrings], dtype=bool)
print(f"Qubit 0 measured |1> in {qubit0_results.sum()} out of {len(qubit0_results)} shots")

# Compute correlator: did qubit 0 and qubit 1 agree?
qubit1_results = np.array([int(bs[-2]) for bs in bitstrings], dtype=bool)
agreement = (qubit0_results == qubit1_results).mean()
print(f"Qubit 0 and 1 agree {agreement:.1%} of the time")
# Expected: ~100% for a GHZ state

Estimator: Expectation Values

The Estimator runs parameterized circuits against observables and returns the expectation value <O> with a standard error estimate.

Parameterized RY Circuit and Z Expectation Value

import numpy as np
from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import StatevectorEstimator

# Parameterized single-qubit circuit
theta = ParameterVector("theta", 1)
qc = QuantumCircuit(1)
qc.ry(theta[0], 0)

# Observable: measure Z on qubit 0
obs = SparsePauliOp("Z")

# Sweep theta from 0 to 2*pi
estimator = StatevectorEstimator()
angles = np.linspace(0, 2 * np.pi, 50)

# Build PUBs: one per parameter value
pubs = [(qc, obs, [[angle]]) for angle in angles]

job = estimator.run(pubs)
result = job.result()

# Extract expectation values and standard errors
exp_values = [result[i].data.evs for i in range(len(angles))]
std_errs   = [result[i].data.stds for i in range(len(angles))]

print(f"At theta=pi/2: <Z> = {exp_values[12][0]:.4f}")
# Expected: <Z> ≈ 0.0  (|+y> state is Z-eigenvalue 0 on average)

SparsePauliOp Complete Reference

SparsePauliOp is the standard way to define observables for the Estimator. It represents a linear combination of Pauli operators.

Single-Qubit Observables

from qiskit.quantum_info import SparsePauliOp

# Single Pauli Z on one qubit
obs_z = SparsePauliOp("Z")

# Single Pauli X
obs_x = SparsePauliOp("X")

# Linear combination: 0.5*Z + 0.3*X
obs_combo = SparsePauliOp.from_list([("Z", 0.5), ("X", 0.3)])
print(obs_combo)
# SparsePauliOp(['Z', 'X'], coeffs=[0.5, 0.3])

Two-Qubit Correlators

# ZZ correlator: measures spin-spin correlation
obs_zz = SparsePauliOp("ZZ")

# XX + YY for transverse coupling
obs_xy = SparsePauliOp.from_list([("XX", 1.0), ("YY", 1.0)])

Heisenberg Hamiltonian on 3 Qubits

The isotropic Heisenberg model on 3 qubits with coupling constant J:

H = J * (X₀X₁ + Y₀Y₁ + Z₀Z₁ + X₁X₂ + Y₁Y₂ + Z₁Z₂)

J = 1.0

# Pauli strings act on all qubits. "IXX" means I on qubit 2, X on qubit 1, X on qubit 0.
# Qiskit uses little-endian ordering: rightmost character is qubit 0.
heisenberg_3q = SparsePauliOp.from_list([
    ("IXX", J), ("IYY", J), ("IZZ", J),   # qubits 0-1 interaction
    ("XXI", J), ("YYI", J), ("ZZI", J),   # qubits 1-2 interaction
])
print(heisenberg_3q)

Energy Offset with Identity

When modeling molecular Hamiltonians, there is often a constant energy offset (nuclear repulsion energy). Include it with the identity operator:

# H = -0.5 * Z0 + 0.3 * Z1 + 0.2 * Z0Z1 + 0.7 * I
mol_hamiltonian = SparsePauliOp.from_list([
    ("IZ", -0.5),   # Z on qubit 0
    ("ZI",  0.3),   # Z on qubit 1
    ("ZZ",  0.2),   # ZZ correlator
    ("II",  0.7),   # identity (energy offset)
])

VQE with Estimator: H2 Ground State Energy

This complete example uses the Estimator to find the ground state energy of molecular hydrogen (H2) using the Variational Quantum Eigensolver (VQE) algorithm.

The qubit Hamiltonian for H2 at equilibrium bond length (0.735 Angstroms), obtained from a Hartree-Fock calculation followed by a parity mapping with two-qubit reduction, is:

import numpy as np
from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import StatevectorEstimator
from scipy.optimize import minimize

# H2 qubit Hamiltonian (parity mapping, 2-qubit reduction)
H2_op = SparsePauliOp.from_list([
    ("II", -1.0523732),
    ("IZ",  0.3979374),
    ("ZI", -0.3979374),
    ("ZZ", -0.0112801),
    ("XX",  0.1809312),
])

# Hardware-efficient ansatz: RY-CNOT-RY on 2 qubits
params = ParameterVector("theta", 4)
ansatz = QuantumCircuit(2)
ansatz.ry(params[0], 0)
ansatz.ry(params[1], 1)
ansatz.cx(0, 1)
ansatz.ry(params[2], 0)
ansatz.ry(params[3], 1)

# Set up the Estimator
estimator = StatevectorEstimator()

# Cost function: returns <H> for given parameters
def cost_fn(theta_vals):
    pub = (ansatz, H2_op, [theta_vals.tolist()])
    job = estimator.run([pub])
    result = job.result()
    energy = result[0].data.evs[0]
    return energy

# Classical optimization using COBYLA
x0 = np.random.uniform(-np.pi, np.pi, size=4)
opt_result = minimize(cost_fn, x0, method="COBYLA", options={"maxiter": 300})

print(f"Optimized energy: {opt_result.fun:.6f} Hartree")
print(f"Exact ground state: -1.857275 Hartree")
print(f"Error: {abs(opt_result.fun - (-1.857275)):.6f} Hartree")
# The optimized energy should converge near -1.857 Hartree

The exact ground state energy of this H2 qubit Hamiltonian is approximately -1.857275 Hartree. This specific set of coefficients comes from an early experimental VQE study using a particular basis set and parity mapping. A well-converged VQE run with this 4-parameter ansatz typically reaches within 0.001 Hartree of the exact value, since the ansatz has enough expressibility to represent the ground state for this 2-qubit system.

Parameterized Circuits and Batching

When you need to evaluate a parameterized circuit at many different parameter values, batching all evaluations into a single .run() call is dramatically faster than submitting them one at a time.

Batching 50 Parameter Sets in One Call

import time
import numpy as np
from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import StatevectorEstimator

# Parameterized ansatz
params = ParameterVector("p", 2)
qc = QuantumCircuit(2)
qc.ry(params[0], 0)
qc.ry(params[1], 1)
qc.cx(0, 1)

obs = SparsePauliOp("ZZ")
estimator = StatevectorEstimator()

# Generate 50 random parameter sets
rng = np.random.default_rng(42)
all_params = rng.uniform(-np.pi, np.pi, size=(50, 2))

# --- Batched: one call with 50 PUBs ---
t0 = time.perf_counter()
pubs = [(qc, obs, [p.tolist()]) for p in all_params]
job = estimator.run(pubs)
result = job.result()
batched_values = [result[i].data.evs[0] for i in range(50)]
t_batched = time.perf_counter() - t0

# --- Sequential: 50 separate calls ---
t0 = time.perf_counter()
sequential_values = []
for p in all_params:
    job = estimator.run([(qc, obs, [p.tolist()])])
    result = job.result()
    sequential_values.append(result[0].data.evs[0])
t_sequential = time.perf_counter() - t0

print(f"Batched:    {t_batched:.3f}s")
print(f"Sequential: {t_sequential:.3f}s")
print(f"Speedup:    {t_sequential / t_batched:.1f}x")

On a local simulator, the speedup is modest (2-5x) because the overhead is just Python function-call latency. On real hardware, the difference is enormous: batched PUBs become a single job in the queue, while sequential calls create 50 separate jobs, each waiting in the queue independently. This can mean the difference between 5 minutes and 5 hours of wall-clock time.

Session-Based Execution for Hardware

A Qiskit Runtime Session is a persistent connection to a specific backend that keeps the QPU reserved for your workload between job submissions. Without a session, each .run() call enters the queue independently. With a session, subsequent jobs skip the queue and execute immediately on the reserved hardware.

Why Sessions Matter

• Reduced queue wait time. Your second, third, and fourth jobs do not re-enter the queue.
• Consistent calibration. All jobs in a session run on the same calibration cycle, so systematic errors are consistent across measurements.
• Interactive workflows. Variational algorithms (VQE, QAOA) submit many short jobs in a loop. Sessions make this practical on shared hardware.

With-Statement Pattern (Recommended)

from qiskit_ibm_runtime import QiskitRuntimeService, Session, SamplerV2, EstimatorV2

service = QiskitRuntimeService()
backend = service.least_busy(operational=True, simulator=False)

with Session(service=service, backend=backend) as session:
    # All primitives created inside the session share the reservation
    sampler = SamplerV2(session=session)
    estimator = EstimatorV2(session=session)

    # Job 1: sample a Bell state
    job1 = sampler.run([(isa_bell_circuit, [], 4096)])
    counts = job1.result()[0].data.meas.get_counts()
    print("Bell counts:", counts)

    # Job 2: estimate an expectation value
    job2 = estimator.run([(isa_circuit, isa_obs, [[0.5]])])
    ev = job2.result()[0].data.evs
    print("Expectation value:", ev)

    # Job 3: another estimation with different parameters
    job3 = estimator.run([(isa_circuit, isa_obs, [[1.0]])])
    ev2 = job3.result()[0].data.evs
    print("Expectation value 2:", ev2)

# Session closes automatically here; QPU reservation is released

Manual Open/Close Pattern

from qiskit_ibm_runtime import QiskitRuntimeService, Session, SamplerV2

service = QiskitRuntimeService()
backend = service.least_busy(operational=True, simulator=False)

session = Session(service=service, backend=backend)

try:
    sampler = SamplerV2(session=session)

    for i in range(10):
        job = sampler.run([(isa_circuit, [], 1024)])
        result = job.result()
        print(f"Iteration {i}: {result[0].data.meas.get_counts()}")

finally:
    session.close()  # Always close to release the QPU

Switching to Real Hardware

The code above runs locally. To target an IBM backend, swap the import and wrap the primitive with the backend:

from qiskit_ibm_runtime import QiskitRuntimeService, EstimatorV2
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

service = QiskitRuntimeService()
backend = service.least_busy(operational=True, simulator=False)

# Transpile the circuit for the target backend
pm = generate_preset_pass_manager(optimization_level=1, backend=backend)
isa_circuit = pm.run(qc)

# Map the observable to the ISA layout
isa_obs = obs.apply_layout(isa_circuit.layout)

estimator = EstimatorV2(backend=backend)
job = estimator.run([(isa_circuit, isa_obs, [[np.pi / 2]])])
result = job.result()
print(result[0].data.evs)

The transpilation step is explicit here because real hardware requires ISA (Instruction Set Architecture) circuits that conform to native gates and qubit connectivity. The apply_layout() call aligns the observable with the physical qubit mapping the transpiler chose.

Batching Multiple PUBs

Both primitives accept a list of PUBs in a single .run() call. This reduces per-job overhead on the hardware queue:

from qiskit import QuantumCircuit
from qiskit.primitives import StatevectorSampler

# Two independent circuits
qc_bell = QuantumCircuit(2)
qc_bell.h(0)
qc_bell.cx(0, 1)
qc_bell.measure_all()

qc_ghz = QuantumCircuit(3)
qc_ghz.h(0)
qc_ghz.cx(0, 1)
qc_ghz.cx(0, 2)
qc_ghz.measure_all()

sampler = StatevectorSampler()

# Submit both in one job; PUB format: (circuit, parameter_values, shots)
job = sampler.run([(qc_bell, [], 2048), (qc_ghz, [], 2048)])
result = job.result()

bell_counts = result[0].data.meas.get_counts()
ghz_counts  = result[1].data.meas.get_counts()
print(bell_counts)  # {'00': ~1024, '11': ~1024}
print(ghz_counts)   # {'000': ~1024, '111': ~1024}

Each PUB result is accessed by index in the returned PrimitiveResult.

Classical Register Naming and Accessing Results

When a circuit has multiple classical registers, the BitArray for each register is accessed by its name. This is common when you measure different groups of qubits for different purposes.

from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister
from qiskit.primitives import StatevectorSampler

qr = QuantumRegister(3, "q")
cr_a = ClassicalRegister(2, "meas_a")
cr_b = ClassicalRegister(1, "meas_b")

qc = QuantumCircuit(qr, cr_a, cr_b)
qc.h(0)
qc.cx(0, 1)
qc.cx(0, 2)

# Measure qubits 0 and 1 into register 'meas_a'
qc.measure(qr[0], cr_a[0])
qc.measure(qr[1], cr_a[1])

# Measure qubit 2 into register 'meas_b'
qc.measure(qr[2], cr_b[0])

sampler = StatevectorSampler()
result = sampler.run([(qc, [], 2048)]).result()

# Access each register by name
counts_a = result[0].data.meas_a.get_counts()
counts_b = result[0].data.meas_b.get_counts()

print("Register meas_a:", counts_a)
# {'00': ~1024, '11': ~1024}

print("Register meas_b:", counts_b)
# {'0': ~1024, '1': ~1024}

The attribute name on result[0].data matches the ClassicalRegister name you defined. If you use measure_all(), the default register name is "meas".

Transpilation Optimization Levels

When targeting real hardware, the pass manager’s optimization_level controls how aggressively the circuit is rewritten before execution:

Level	What it does
0	No optimization. Maps to hardware trivially. Fastest compile time.
1	Light optimization: cancels redundant gates, trivial 1Q reductions. Default.
2	Heavier optimization: commutation analysis, unitary synthesis.
3	Heaviest: full unitary synthesis and noise-aware layout/routing.

For NISQ research, level 1 or 2 is most common. Level 3 can significantly reduce two-qubit gate count at the cost of longer compilation.

Error Mitigation with Estimator

EstimatorV2 supports built-in error mitigation through the resilience_level option. Each level enables progressively more sophisticated mitigation at the cost of additional circuit executions.

Resilience Level 0: No Mitigation

Raw hardware output with no post-processing. Use this when you need the fastest possible execution or when you are benchmarking hardware noise directly.

from qiskit_ibm_runtime import EstimatorV2, EstimatorOptions

options = EstimatorOptions()
options.resilience_level = 0

estimator = EstimatorV2(backend=backend, options=options)
job = estimator.run([(isa_circuit, isa_obs, [[np.pi / 4]])])

Overhead: None. One circuit execution per PUB.

Resilience Level 1: Twirling + Readout Mitigation

Applies Pauli twirling to convert coherent errors into stochastic noise (which is easier to model and correct) and applies readout error mitigation to correct measurement errors.

options = EstimatorOptions()
options.resilience_level = 1

estimator = EstimatorV2(backend=backend, options=options)
job = estimator.run([(isa_circuit, isa_obs, [[np.pi / 4]])])

Overhead: Low. Requires a small number of additional calibration circuits for readout error characterization, but the total circuit count is only modestly higher than level 0.

When to use: This is the best default for most experiments. The accuracy improvement is substantial relative to the cost.

Resilience Level 2: Zero Noise Extrapolation (ZNE)

Runs the circuit at several artificially amplified noise levels, then extrapolates to the zero-noise limit. This provides the most accurate expectation values available on current hardware.

options = EstimatorOptions()
options.resilience_level = 2

estimator = EstimatorV2(backend=backend, options=options)
job = estimator.run([(isa_circuit, isa_obs, [[np.pi / 4]])])

Overhead: Significant. ZNE typically requires 3-5x more circuit executions because each PUB is run at multiple noise amplification factors. For example, with three noise levels, the circuit count triples.

When to use: When you need the most accurate expectation values possible, such as for final VQE energy estimates or precise Hamiltonian simulation results. Avoid using it during the optimization loop of a variational algorithm (use level 1 instead, then do a final evaluation at level 2).

Pauli Twirling Configuration

Pauli twirling randomizes the noise channel of two-qubit gates, converting coherent errors into incoherent (stochastic) Pauli errors. You can configure twirling independently of the resilience level:

options = EstimatorOptions()
options.resilience_level = 1
options.twirling.enable_gates = True    # twirl two-qubit gates
options.twirling.enable_measure = True  # twirl measurements
options.twirling.num_randomizations = 32  # number of random twirling instances

estimator = EstimatorV2(backend=backend, options=options)

More randomizations produce a cleaner conversion of coherent noise to stochastic noise, but each randomization adds circuit executions. For most use cases, 32 randomizations provides a good balance.

Comparing StatevectorSampler vs. Noisy Simulation

Both StatevectorSampler and the Aer simulator run locally, but they serve different purposes.

StatevectorSampler performs exact statevector simulation. There is no noise, no decoherence, and no gate errors. This is ideal for verifying circuit correctness and developing algorithms.

AerSimulator supports noise models that replicate the error characteristics of real hardware. This lets you estimate how your circuit will perform on a specific backend before consuming real QPU time.

Adding a Fake Backend Noise Model

from qiskit import QuantumCircuit
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime.fake_provider import FakeSherbrooke

# Create a noisy simulator that mimics the Sherbrooke backend
fake_backend = FakeSherbrooke()
noisy_sim = AerSimulator.from_backend(fake_backend)

# Build and transpile a circuit for this fake backend
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

qc = QuantumCircuit(2)
qc.h(0)
qc.cx(0, 1)
qc.measure_all()

pm = generate_preset_pass_manager(optimization_level=1, backend=fake_backend)
isa_circuit = pm.run(qc)

# Run with the noisy simulator using SamplerV2
from qiskit_ibm_runtime import SamplerV2

sampler = SamplerV2(backend=noisy_sim)
job = sampler.run([(isa_circuit, [], 4096)])
result = job.result()

counts = result[0].data.meas.get_counts()
print(counts)
# You will see small counts for '01' and '10' due to noise,
# unlike the exact StatevectorSampler which only produces '00' and '11'

This workflow is valuable for noise-aware algorithm development. You can iterate locally with realistic noise before submitting to the actual QPU.

Migration Guide: execute() and backend.run() to Primitives

If you have existing code that uses the deprecated execute() or backend.run() API, here is how to migrate to primitives.

Old Pattern (Deprecated)

# OLD WAY (deprecated, will be removed)
from qiskit import QuantumCircuit, transpile

qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

# backend = ...  (some backend instance)
transpiled = transpile(qc, backend)
job = backend.run(transpiled, shots=1024)
result = job.result()
counts = result.get_counts()    # returns a plain dict
print(counts)
# {'00': 512, '11': 512}

New Pattern (Primitives)

# NEW WAY (current API)
from qiskit import QuantumCircuit
from qiskit.primitives import StatevectorSampler

qc = QuantumCircuit(2, 2)
qc.h(0)
qc.cx(0, 1)
qc.measure([0, 1], [0, 1])

sampler = StatevectorSampler()
job = sampler.run([(qc, [], 1024)])
result = job.result()
counts = result[0].data.c.get_counts()   # BitArray, not a plain dict
print(counts)
# {'00': ~512, '11': ~512}

Migration Checklist

# Migration checklist:
# 1. Replace `from qiskit import execute` with primitive imports
# 2. Replace `backend.run(circuit)` with `sampler.run([(circuit,)])`
# 3. Replace `result.get_counts()` with `result[0].data.<register>.get_counts()`
# 4. If you need expectation values, switch from Sampler to Estimator
# 5. Replace manual transpile() calls with generate_preset_pass_manager()
# 6. For observables, call obs.apply_layout(isa_circuit.layout) after transpilation
# 7. Wrap hardware primitives in a Session for multi-job workflows
# 8. Set resilience_level in EstimatorOptions for error mitigation

Common Mistakes

1. Forgetting to Transpile Before Using EstimatorV2 on Real Hardware

The runtime primitives on real hardware expect ISA circuits (circuits composed of the backend’s native gates on its physical qubits). If you pass an untranspiled circuit, you get an error.

# WRONG: passing an abstract circuit to hardware
estimator = EstimatorV2(backend=backend)
job = estimator.run([(qc, obs, [[0.5]])])  # Error: circuit not in ISA form

# CORRECT: transpile first
pm = generate_preset_pass_manager(optimization_level=1, backend=backend)
isa_circuit = pm.run(qc)
isa_obs = obs.apply_layout(isa_circuit.layout)
job = estimator.run([(isa_circuit, isa_obs, [[0.5]])])

2. Confusing BitArray with dict

If you write counts = result.get_counts() like the old API, you get an AttributeError. The result object is a PrimitiveResult, and you must index into it.

# WRONG
counts = result.get_counts()

# CORRECT
counts = result[0].data.meas.get_counts()

3. Wrong PUB Format (Missing List Wrapping)

The .run() method expects a list of PUBs, even if you only have one. Passing a bare tuple causes confusing errors.

# WRONG: bare tuple, not wrapped in a list
job = sampler.run((qc, [], 1024))

# CORRECT: list of PUBs
job = sampler.run([(qc, [], 1024)])

4. Not Calling apply_layout() on Observables After Transpilation

After transpilation, the logical-to-physical qubit mapping may change. If you pass the original observable without adjusting for the new layout, the Estimator measures the wrong qubits.

# WRONG: observable still references logical qubits
pm = generate_preset_pass_manager(optimization_level=1, backend=backend)
isa_circuit = pm.run(qc)
job = estimator.run([(isa_circuit, obs, [[0.5]])])  # obs is on wrong qubits

# CORRECT: apply the layout
isa_obs = obs.apply_layout(isa_circuit.layout)
job = estimator.run([(isa_circuit, isa_obs, [[0.5]])])

5. Using Sampler When You Need Estimator (and Vice Versa)

The Sampler returns bitstring distributions. The Estimator returns expectation values. Choosing the wrong one creates unnecessary work.

• Use the Sampler when you need to see which bitstrings the circuit produces (quantum state tomography, combinatorial optimization readout, random circuit sampling).
• Use the Estimator when you need the numerical value of an observable (energy estimation in VQE, correlation functions, order parameters in quantum simulation).

A common antipattern is using the Sampler to compute <Z> by manually counting bitstrings and averaging. The Estimator does this automatically, with built-in error mitigation, and returns the standard error.

6. Closing a Session Prematurely in a Loop

If you open and close a session inside a loop, each iteration creates a new session (and a new queue entry), defeating the purpose.

# WRONG: new session every iteration
for params in param_list:
    with Session(service=service, backend=backend) as session:
        estimator = EstimatorV2(session=session)
        job = estimator.run([(isa_circuit, isa_obs, [params])])
        print(job.result()[0].data.evs)

# CORRECT: one session wrapping the entire loop
with Session(service=service, backend=backend) as session:
    estimator = EstimatorV2(session=session)
    for params in param_list:
        job = estimator.run([(isa_circuit, isa_obs, [params])])
        print(job.result()[0].data.evs)

Summary

Primitive	Input	Output	Use for
Sampler	circuits + shots	BitArray (bitstring distributions)	Sampling, circuit output analysis
Estimator	circuits + observables + params	Expectation values + std errors	VQE, QAOA cost, Hamiltonian simulation

The primitives API is now the standard interface for IBM Quantum. Any code written against StatevectorSampler or StatevectorEstimator runs identically on real hardware by swapping to SamplerV2 or EstimatorV2 with a backend argument, making it straightforward to prototype locally and deploy to hardware.