Interpretation of occupation numbers and creation-annihilation operators in strongly correlated systems I hope I am not asking anything stupid, but I am having a hard time interpreting some
(seemingly simple) results.
The simplest approximate form of an $N$-electron wave function (w.f.) is a Slater-determinant:
$$ |\Phi\rangle=\prod_{i=1}^{N}a_i^{\dagger}|0\rangle \ , $$
where $|0\rangle$ is the zero electron vacuum state and $a_i^{\dagger}$ is the creation
operator of orbital $\phi_i$.
The corresponding 1-particle reduced density matrix (1RDM) is diagonal:
$$ \gamma_{pq}=\langle a_{q}^{\dagger}a_{p}\rangle=\delta_{pq}n_{p} \ . $$
The $n_p$ occupation numbers are either one for occupied orbitals, or zero for unoccupied 
(virtual) orbitals.
The annihilation operator $a_k$ has the interpretation that it destroys an electron (or 
equivalently, creates a hole), thus mapping the $N$-electron w.f. to one with $(N-1)$ 
electrons.
So, based on the integer character of occupation numbers (say $n_k=1$), we can say that 
every electron is represented by a single orbital ($\phi_k$), and they can be destroyed by the
action of the corresponding annihilation operator ($a_k$).
So far, so good.
However, in a strongly correlated system (e.g. molecules with stretched bonds, or transition metal compounds) the above single-determinantal approximation is known to break down, and one needs a 
better w.f. form, which is usually parametrized as a linear combination of Slater-determinants:
$$ |\Psi\rangle=\sum_{K \in \Omega}C_{K}|\Phi_K\rangle \ .$$
Here $\Omega$ is some finite subspace of the basis spanned by all possible determinants.
Now, $|\Psi\rangle$ is much more suitable for the treatment of strong correlation, but it is 
also harder to interpret (at least for me).
Working on the basis of the so-called natural orbitals, the diagonal form of $\gamma_{pq}$ can
be retained, but at the cost of having fractional occupation numbers with $n_p \in [0,1]$.
My interpretation for this was the following: unlike in the single-determinantal picture, now
an electron is represented by more than one orbital (e.g. if $n_p+n_q+n_r=1$, than orbitals
$\phi_p$, $\phi_q$ and $\phi_r$ represent the electron together). 
In this way $n_p$ can be thought of as an "average occupation", or equivalently, it measures
how much $\phi_p$ contributes to the occupation of the electron.
Note that I am not speaking about fractional electrons physically, it is just an
interpretation.
Now let us act on $|\Psi\rangle$ with $a_p$, and then read off the
number of electrons with the number operator $\hat{N}=\sum_{p}a_{p}^{\dagger}a_p$.
By using $[\hat{N},a_p]=-a_p$, we get
$$ \hat{N}a_p|\Psi\rangle=(N-1)a_p|\Psi\rangle \ .  $$
Alternatively, one could get this result by writing $|\Psi\rangle$ as
$$ |\Psi\rangle=\sum_{\substack{K \in \Omega \\ \phi_p \in    \Phi_K}}C_{K}|\Phi_K\rangle
               +\sum_{\substack{K \in \Omega \\ \phi_p \notin \Phi_K}}C_{K}|\Phi_K\rangle \ , $$
and noticing that $a_p$ destroys everything in the second term.
Here comes my problem: based on the occupation numbers, $\phi_p$ acts like a fraction of an electron, but $a_p$ still annihilates a complete electron!
Does somebody know a way to reconcile these two (seemingly contradicting) observations?
Or more generally, is there a consistent interpretation of fractional occupation numbers
and creation-annihilation operators in the regime of strong correlation?
Note that I do not think there is anything wrong with the theory.
I just think that I have not found the proper interpretation of these results.
A final remark: I know that fractional occupation numbers arise in many other areas of 
physics/chemistry (e.g. in Density Functional Theory).
However, I am currently not particularly interested in these other aspects (unless if they
help answer the current questions, of course).
 A: You are mostly correct in your analysis and interpretation. I would say that the strongly-correlated eigenstate is a complicated superposition of many free-Hamiltonian eigenstates $|\Phi\rangle$, each of them (in your example) with the same total number of electrons $N$ (this doesn't have to be the case, but this isn't the matter right now). Therefore, when we calculate the expectation value of some occupancy of an orbital of the free-theory, we will get a fractional number as some of the participating wave functions have it occupied and some don't. This is nothing to be surprised by - the free-theory orbitals belong to the free-theory, and while we find it convenient to think about the system using their terms, the system itself doesn't regard them as anything special. It is interacting and strongly correlated, and doesn't have to adhere to our notions of what is an electronic  orbital.
When we act with an annihilation operator $a_p$ on the system, however, we do something quite "violent" to it. Not only do we remove an electron from the orbital $a_p$ for all the participating wave-functions where there was one there, but we also completely annihilate the entire part of the superposition that had an occupancy zero at the $p$-orbital to begin with! The resulting wave function now has a definite occupancy of the $p$-orbital and it is zero. Since it was composed of a superposition of wave-functions with $N$ electrons, which now were either annihilated to zero (if they had $p$-orbital empty) or had one electron removed from them (if they had $p$-orbital full), the resulting total number of electrons of those states that survived must be $N-1$, and $\langle n_p \rangle =0$ (no matter if before that it was fractional).
Note, however, that as the resulting state is not an eigenstate of the fully interacting theory (in general), the time-evolution of the state will make it evolve into a state where $\langle n_p \rangle$ is again fractional, but all the states will still have $N-1$ particles in them.
