There should be an easy way to compute expectation values for operators

[nonhermitian]

What should we add?

Currently there is no way good way to compute expectation values in Qiskit. The StateVector and DensityMatrix objects in quantum_info have expectation_value methods, but there is no general way to do this.

In general expectation values are defined with respect to the operator in question and not the state. To quote Sakurai:

We define the expectation value of $A$ taken with respect to $|a\rangle$

So it makes sense to have the operator as the primary object, and having a method by which one can pass Counts, QuasiDistribution or the like and obtain the desired value. Obtaining these values requires the bit-string representation to be available via those objects efficiently.

There is no unique way to define operators in Qiskit, eg there is are Pauli and PauliOp. So not even sure what flavor of operator is primary to Qiskit as well.

[wyqian1027]

In my opinion, the opflow module is quite useful. You define your operator and then sampler of any QuantumInstance (statevector, QASM, noise, etc), then you get your measurement results.

[jlapeyre]

It's certainly true that referring to the expectation value of an observable (often without reference to the state) is common. But, I don't think that it's the only useful way to discuss the average value of the measurement of a particular observable on a system in a particular state.

In their book, Nielsen and Chuang do not use the term expectation value. They use "average" or "expected value" when discussing both classical and quantum systems as well as probability theory. They do introduce the notation $\langle A \rangle$. Doing a quick search, I find only one place (pg 116) where they use this notation. In fact, they go out of their way to say "the average value for these observables, written in the quantum mechanical $\langle \cdot \rangle$ notation, are:". And in this case it's clear that this is because it makes the expression for the CHSH inequality more clear.

In Qiskit, we do frequently use expectation value. Regardless of the term we use for it, we have to decide whether it belongs to the state or the operator. But that's because Sapir-Whorf applies here. Our world view is so constrained by class-based OO languages that we tend to project the constraints from the programming language onto the subject domain.

But there is a way out! There are several multiple-dispatch libraries available for Python. Multiple-dispatch is well suited to mathematical and scientific domains. The problem at hand is a good example.

We would define a function expectation that is not a method that belongs to any class. Several methods would be defined for expectation (in general in different files) like this

@dispatch(PauliOp, StateVector)
def expectation(op, state):
...

@dispatch(PauliOp, DensityMatrix)
def expectation(op, state):
...

@dispatch(PauliOp, QuasiDistribution)
def expectation(op, state):
...

@dispatch(Operator, DensityMatrix)
def expectation(op, state):
...

The consumer would use them like this:

from qiskit.quantum_info import expectation

op = ...
state = ...
expectation(op, state)

[nonhermitian]

@jlapeyre Hmm, I think expectation value and average mean the same thing here. I guess your advocating for less physics orientated terminology? Even the Mike and Ike quote above does talk about the "average value" as corresponding to the operators and not of the states. This is to be expected as the operator is the thing that tells you what you are measuring (assuming a Hermitian operator) and not the state.

I am not sure Sapir-Whorf is really applicable as mathematics is the universal language we should be using. I think it is more a question of how to accomplish the exact same task. I am advocating we follow the mathematical way, in particular notation wise where things are written $\langle O \rangle$, $\mathbb{E}[O]$, etc, and treat the operator as the primary object and not the state. Although I to have fallen into the state first notation, eg. see #7822.

The state first approach is easy to fall into because Qiskit is built around circuits. Now circuits are just evolution operators, but in practice we often assume that they start from the ground state and thus equate them with producing statevectors. So states, or finite sampling distributions obtained from those states, are primary in the eyes of much of Qiskit as well. So perhaps states is the way to go, and what I am advocating for is something that exists on a level above Qiskit Terra and which calls expectation values using states as the primary object. Regardless, this is a big gap in the Qiskit core functionality, and is something that should be added regardless of the implementation details

[jlapeyre]

@nonhermitian, I suppose I do favor the less physics-y terminology for QI and QC. Sakurai and other texts of that period made a conscious break with earlier presentations, such as omitting a long section on the development of QM. I think Mike and Ike did a good job of adapting the presentation of the material to a new context. But, I'm not advocating purifying the Qiskit documentation to use only one approach. Mostly, I'm arguing that none of these different approaches to presentation and notation are right or wrong in the sense that math can be right or wrong.

In some physics writing, and more so in math, you find things like $\langle O \rangle_\phi$ and $\mathbb{E}_p[O]$ as the fundamental notation, with the state or measure omitted if it is clear from the context. (For example the wikipedia article uses this notation.)

Since I'm appealing to authorities, I'll mention that Mark Wilde's book has

It is common for physicists to denote the expectation as $\langle Z \rangle \equiv \langle \psi | Z | \psi \rangle$ when it is understood that the expectation is with respect to the state $|\psi\rangle$.

In any case, I'd argue that a reasonable way to use math notation to present the expected value (for pure states) with the standard physics notation (provided $\langle \phi | A | \psi \rangle$ is already defined) is:

Let $X$ be a complex Euclidean space and ${\cal H}(X)$ the Hermitian operators on this space. The expected value is a function $$\langle \cdot \rangle_\cdot : {\cal H}(X) \times X \rightarrow \mathbb{R}; (O, |\phi\rangle) \mapsto \langle O \rangle_\phi $$ given by $\langle O \rangle_\phi = \langle \phi | O | \phi \rangle$.

The main point is that the operator, the state, and the resulting value enter on equal footing, just as with any other map that maps two values to one value.

Mentioning Sapir-Whorf is probably an overstatement; I just meant that a class-based-OO (CBOO) approach is so dominant that we don't imagine other options. Or at least we'd say "Since we have a big software project that is already committed to a CBOO approach, we have to decide whether expectation value belongs to our operator classes or our state classes." Of course, you could do both, but that would be confusing.

In any case, you have been bringing up practical needs in several issues. If one solution ruins our abstractions, we should find another. As I said above, I think making mathematical operations belong to one object in an expression is forcing CBOO thinking on a problem that does not naturally have this structure. It doesn't really naturally have a multiple-dispatch structure either. But, I think MD is a much better fit. It does add another dependency, but it is pure python and quite light. I'm not 100% sure about how perfect it is for this situation, I'd have to look at details. But, I'd be willing to make a trial PR. (I made a PR in the past, putting MD in several places, just to show what it looks like.)

[wyqian1027]

The consumer would use them like this:

from qiskit.quantum_info import expectation

op = ...
state = ...
expectation(op, state)

@jlapeyre I think I understand what you mean now and this is also what I prefer to think, expectation as a method that works together with state and operator. In practice, I usually go by using CircuitSampler to convert an expectation operator to expectation values, where the expectation operator is built by composing the initial state and the ansatz. So, in a way, the workflow is similar to what you described above.

From a physics perspective, the expectation value, $\braket{\psi|\hat{O}|\psi}$, is always on the quantum state itself, i.e., $\ket{\psi}$, be it omitted or not. So I think it is fine that functions like "Statevector" and "DensityMatrix" have those convenient functions that take an operator to compute expectation values, which is helpful for directly manipulating your simulation results. Though from a software design point of view, it might be a different story.

[nonhermitian]

I do not think the operator and the state are on the same footing exactly. Namely it is the operator that determines much of the physical properties of the expectation value. The most obvious of course is the Hermitian structure, or lack thereof, determines whether the expectation value corresponds to a physical observable or not. The state has no such bearing on that.

From a programmatic viewpoint, the operator instance makes more sense as the place to put a method. For example, Qiskit has Pauli, PauliOp , PauliSumOp etc etc, and each would have its expectation value computed differently. If this was done in the distributions then it would have to know about all these types, and dispatch appropriately. In contrast, each operator type can carry around its own definition in its expval method.

[nonhermitian]

I am also not a big fan of the expectation(op, state) way of doing things as it is a very Matlab way of doing things and does not leverage the OO nature of Python.

[ikkoham]

I believe that Estimator in Primitives provides the unified interface. But we can have another approach for quantum_info.

First, Operator class cannot have expval methods because it introduces cyclic dependency. (Many deprecations of methods in states are necessary...) States and observables are the dual concept. There are no priority. so I vote expectation function John proposed. (Either MD would be fine. I think you would say it isn't clean, but we could branch out with isinstance in the fuction.)

[nonhermitian]

The Estimator doesn't actually provide a unified interface. If all someone wanted was expectation values then yes. But what if I want to compute expectation values and also plot the distribution that generated that value, or compute the fidelity with the distribution? Then I need to drop down to the Sampler and then I am back to needing an expectation value routine.

Operator class cannot have expval methods because it introduces cyclic dependency

This I do not understand. Given an operator and a collection of bit-strings, nothing but iterating over the entries in the collection is needed to compute the expectation value. Perhaps your saying the refactoring needed to do so is difficult?

[ikkoham]

I thought the purpose was to calculate the expectation value in this issue (but I found it was wrong...). Also, in discussion there are two types of expectation values. (we are confusing.) One is the quantum expectation $<\psi |A \psi>$ that John and I are talking about, and the other is the Z-diagonal expectation (in a sense, the classical expectation) that is calculated from the probability distribution you mention. If you are arguing that a method to calculate the classical expectation is needed for Counts or Probability, I completely agree.

The point you don't understand comes from the difference in expectations as described above, and I mistakenly thought we were talking about quantum expectations. My apologies.

[ikkoham]

We are considering adding a DiagonalOperator class to quantum_info for the next release. For classical expectation, I think it is sufficient to support SparsePauliOp and DiagonalOperator. Other operator classes can be easily converted to them.

I would like to propose following methods:

SparsePauliOp.classical_expectation(data :Counts | Probability | QuasiDistribution)
DiagonalOperator.classicall_expectation(data :Counts | Probability | QuasiDistribution)

[jlapeyre]

@ikkoham , I'm not sure that @nonhermitian was only concerned with computing expectations of things like Counts, although that was emphasized as a pain point. I was reading this as arguing first that the operators are the natural place to put a method for expected value. And an important case is when the state is Counts, etc.

I don't see a large distinction between what you call classical and quantum expectations. If you assign a number $-1^{n_1}$ to each bitstring where $n_1$ is the number of 1s in a bitstring, and average them, you have computed the value of $Z\otimes\ldots\otimes Z$ on the state from which the counts were collected, right? We may want to organize this in software differently than other expectation values. But, it's still the expectation value of a quantum state.

But, people have also mentioned in other issues (https://github.com/Qiskit/qiskit-terra/issues/7822) the following:

• Counts are not documented as being the result of a Z measurement. In principle a backend could support other bases. In fact, we could do the same in software. Provide an abstraction saying that we want to measure in the $X$ basis, then insert a layer of $H$, then collect counts. And we could tag the the resulting data structure by the basis it was measured in.
• Since the phase information is no longer in the counts (or quasiprobability), you can only have I and Z in your measurement string.

I suppose we could offer a method that computes the expectation with $IAIBCII$, etc. where $A,B,C\ldots$ are not specified, but the user understands that they correspond to whatever measurement basis was used. The input would just be the positions without a I. If I understand correctly, this is the most general expectation you can compute from the probabilities and the implementation doesn't depend on the measurement basis. But the semantics of method are that it computes the (partial) expectation value in whatever basis you measured in. It seems that a method for this operation must belong to Counts and not to the observable.

When to use MD: If it is just a few if isinstance here and there, it's not worth the trouble. But, I think there is the possibility to find a situation where it really removes some complexity. This may or may not meet the bar. In particular, on a PL level, the semantics of the expectation value are rather different between: 1) That outlined above (for Counts) , and 2) that of the expectation given an explicit numeric matrix and vector representing an observable and state. Enough so, that trying to fit them into the same generic function might be stretching the meaning of that function too much.

EDIT: My proposal is then a method

Counts.expectation(id_positions : List)

This computes the expectation value in whatever basis the measurement was taken on whatever state the system was in. We substitute I at all bit positions given in id_positions. And id_position = [] means compute the expectation of the full measurement operator.

EDIT 2: An alternative is that id_positions be the complement of the definition give above. Then counts.expectation([3,6]) would be the expectation value of $A_3 A_6$, where $A_i$ is the measurement basis used on the $i$ th qubit.

[jyu00]

The Estimator doesn't actually provide a unified interface. If all someone wanted was expectation values then yes. But what if I want to compute expectation values and also plot the distribution that generated that value, or compute the fidelity with the distribution? Then I need to drop down to the Sampler and then I am back to needing an expectation value routine.

Estimator was specifically created to return expectation values. If the concern is it doesn't give you the distribution, can't we just add that to the estimator result?

[jlapeyre]

If I understand what's asked for here, I think Estimator is probably not the right tool. The documentation and implementation require that the observables be represented by lists of circuits and those circuits are required to have parameters. (They don't really require parameters, but in this case you are required to pass an empty list of values when calling). Estimator looks convenient for wrapping a lot of workflows that come up in algorithms.

I'm thinking of a simple, uniform, interface that takes something representing an observable and something representing a state and returns the expectation value.

[nonhermitian]

I would also say that Qiskit is about more than the Runtime. Other hardware vendors might not implement the Runtime interfaces, and anyway the base implementations are statevector only. Thus one still needs the expectation values from counts functionality somewhere if the reference implementations are modified to work with counts.

The Estimator could return counts, but then one could ask why do I need the Sampler? You also get into situations where I run readout mitigation, but cannot return the quasi-distribution. So one has to return the original. But then for parity reasons one should probably return the original distribution in the Sampler as well. And then people will complain about the data sizes being too large, etc etc. So taken as a whole, it seems like a lot of round-a-bout to avoid making a method whose functionality we already know how to code up, but have not done so to date.

[jyu00]

My point was can't we improve qiskit-terra's implementation of Estimator to be the "simple, uniform interface" that we want? Qiskit is about more than Runtime, but the primitive interfaces live in qiskit-terra. And the primitives are advertised as ways to perform foundational tasks such as calculating expectation values. If there is a presumably better method for doing that, then why do we need Estimator in terra? And doesn't this mean I'll need to learn a new interface when trying to do the same thing on ibm or other third party hardware?

[nonhermitian]

So yeah when improving the Estimator your going to end up writing something like counts.expectation_value and then the request of this Issue is fulfilled, and others outside of Runtime can use it as well. But perhaps expanding the Runtime base implementations is the motivation for doing so.

[jlapeyre]

My point was can't we improve qiskit-terra's implementation of Estimator to be the "simple, uniform interface" that we want?

I don't think so. Not in the sense that I meant. By itself, the phrase I used: "simple, uniform interface" is not precise enough to understand what I meant. Here is a bit more detail...

I was not involved in the design of the primitives. But, looking at Estimator, it seems to be a high-level, generic tool for handling complexity in implementing high-level (not just one circuit) algorithms. In particular, hybrid quantum-classical algorithms. To me, it looks like a good attempt at "finding the right abstraction". It has a clear purpose that is (probably) very widely applicable. And the implementation (at least here, with Statevector) is simple. We'll see how this plays out when using it in several algorithms.

But, what I'm thinking about here (and I think it overlaps substantially with what @nonhermitian is talking about) is a lower-level abstraction. Perhaps something that would be used inside implementations of primitives.

It is possible to modify Estimator to do low-level expectation calculations. But, this would introduce complexity at the conceptual as well as the coding level: It would muddy the intended abstraction, and add significantly more branchy code.

[jlapeyre]

I pushed a branch implementing (some) methods for expectation_value using multiple dispatch.

Following is an example script that calls these methods:

import numpy as np

from qiskit.quantum_info.states import Statevector
from qiskit.quantum_info.operators import Pauli
from qiskit.quantum_info import expectation_value

diagop = Pauli('ZZ')
oneop = Pauli('Z')

state = Statevector([1, 2, 2, 1])
state = state / np.linalg.norm(state.data)


cs = state.sample_counts(10000)

# Compute expectation of ZZ operator on state
print(expectation_value(diagop, state))

# Explicitly specify all qubits
print(expectation_value(diagop, state, range(len(diagop))))

# Explicitly specify only first qubit
print(expectation_value(oneop, state, [0]))

# Compute expectation value from counts, using all qubits
print(expectation_value(cs))

# Explicitly specify all qubits
# Note, I probably will use a different interface here. Best is probably
# following exactly how qargs works above.
print(expectation_value(cs, '11'))

# Explicitly specify only first qubit
print(expectation_value(cs, '10'))

The meaning of method expectation_value(cs), is found in comments

Compute the expectation value of Counts in the same basis that was measured to collect the counts, for a subset of the qubits. op is a binary string specifying which qubits to include. For example 1100 specifies computing the expectation value of the first two qubits in a four-qubit register.

In addition to aiding in designing a uniform API, the complexity and repetition in the code can be reduced with MD. One example of the benefit of using multiple dispatch in this mathematical setting is the following: This branch implements expectation_value for both Statevector and DensityMatrix in a single method, rather than two. Furthermore, all details of the structure of the operator are hidden in this method. Note that some, but not all, of this can be done with single dispatch; for example with class methods. In any case, I find that the design style with MD encourages the factoring of code shown above.

[nonhermitian]

This is done