Can we write the density operator as a sum of mixed states?

Question

In every resource I find (like Nielsen and Chuang or online courses), the density operator is defined as follows: we consider a sequence of pure states $\left|\psi_i\right\rangle$ with associated probabilities $p_i$. If the system is in state $\left|\psi_i\right\rangle$ with probability $p_i$, then its density operator is given by: $$\rho = \sum_ip_i\left|\psi_i\right\rangle\left\langle\psi_i\right|.$$ However, I do not understand the necessity of $\left|\psi_i\right\rangle$ being pure states for this. For instance, let us say that the system is in state (potentially mixed) $\rho_i$ with probability $p_i$ (that is, $\rho_i$ is a density operator). Then we can prove that its density matrix is given by: $$\rho=\sum_ip_i\rho_i.$$ For this, I use the Lemma defined in this answer, assuming it is correct. We consider an arbitrary unitary $\mathbf{U}$ and an arbitrary basis state $|x\rangle$. Applying $\mathbf{U}$ on the system, we obtain, with probability $p_i$: $$\mathbf{U}\rho_i\mathbf{U}^\dagger.$$ The probability of measuring $|x\rangle$ in this situation is thus given by: $$\mathrm{tr}\left(|x\rangle\langle x|\mathbf{U}\rho_i\mathbf{U}^\dagger\right).$$ Since this situation happens with probability $p_i$, the probability of measuring $|x\rangle$ is finally given by: $$\sum_ip_i\mathrm{tr}\left(|x\rangle\langle x|\mathbf{U}\rho_i\mathbf{U}^\dagger\right).$$ On the other hand, applying $\mathbf{U}$ on $\rho$ gives: $$\mathbf{U}\rho\mathbf{U}^\dagger=\sum_ip_i\mathbf{U}\rho_i\mathbf{U}^\dagger$$ which means that the probability of measuring $|x\rangle$ is given by, using the linearity property of the trace: $$\sum_ip_i\mathrm{tr}\left(|x\rangle\langle x|\mathbf{U}\rho_i\mathbf{U}^\dagger\right).$$ Hence, for any unitary $\mathbf{U}$ and any basis state $|x\rangle$, the probabilities of measuring $|x\rangle$ after having applied $\mathbf{U}$ on those systems are equal, hence their density matrices are identical.

I cannot see the error here, but on the other hand, I find it surprising that I did not find any mention to this (or I've missed it, which is also quite likely). Is this result correct? Otherwise, where's my mistake?

Answer 1

All your working is correct. Of course, if you want to define a density operator/matrix, you cannot define it in terms of another density operator/matrix which is why you will not find an initial presentation of density matrices going in this way. That said, you may want to look up the definition of a mixed separable state. This is typically done in the form $$ \rho_{\text{sep}}=\sum_ip_i\sigma^A_i\otimes \sigma^B_i $$ where $\sigma^A_i$ and $\sigma^B_i$ are density matrices. You could write them just as pure states but you will often see it in this form, consistent with the structure you're talking about.

Answer 2

Due to $\rho$ being a positive and adjoint, you can always spectrally decompose it into a convex combination of pure states, it's eigenvectors. However, since the measurement statistics, or the action of a unitary operator, as you observed, will not change if you do not represent it diagonally, a lot of papers or books will choose not to represent them as this decomposition.