A mixed state is mathematically represented by a bounded, positive trace-class operator with unit trace $\rho : \cal H \to \cal H$.
Here $\cal H$ denotes the complex Hilbert space of the system (it may be nonseparable). The set of mixed states $S(\cal H)$ is a convex body in the complex linear space of trace class operators $B_1(\cal H)$ which is a two-side $*$-ideal of the $C^*$-algebra of bounded operators $B(\cal H)$.
Convex means that if $\rho_1,\rho_2 \in S(\cal H)$ then a convex combination of them, i.e. $p\rho_1 + q\rho_2$ if $p,q\in [0,1]$ with $p+q=1$, satisfies $p\rho_1 + q\rho_2 \in S(\cal H)$.
two-side $*$-ideal means that linear combinations of elements of $B_1(\cal H)$
belong to that space (the set is a subspace), the adjoint of an element of $B_1(\cal H)$ stays in that space as well and $AB, BA \in B_1(\cal H)$ if $A\in B_1(\cal H)$ and $B \in B(\cal H)$.
I stress that, instead, the subset of states $S(\cal H)\subset B_1(\cal H)$ is not a vector space since only convex combinations are allowed therein.
The extremal elements of $S(\cal H)$, namely the elements which cannot be decomposed as a nontrivial convex combinations of other elements, are all of the pure states. They are of the form $|\psi \rangle \langle \psi|$ for some unit vector of $\cal H$. (Notice that, since phases are physically irrelevant the operators $|\psi \rangle \langle \psi|$ biunivocally determine the pure states, i.e. $|\psi\rangle$ up to a phase.)
The space $B_1(\cal H)$ and thus the set $S(\cal H)$ admits at least three relevant normed topologies induced by corresponding norms. One is the standard operator norm $||T||= \sup_{||x||=1}||Tx||$ and the remaining ones are:
$$||T||_1 = || \sqrt{T^*T} ||\qquad \mbox{the trace norm}$$
$$||T||_2 = \sqrt{||T^*T||} \qquad \mbox{the Hilbert-Schmidt norm}\:.$$
It is possible to prove that:
$$||T|| \leq ||T||_2 \leq ||T||_1 \quad \mbox{if $T\in B_1(\cal H)$.}$$
Moreover, it turns out that $B_1(\cal H)$ is a Banach space with respect to $||\cdot||_1$ (it is not closed with respect the other two topologies, in particular, the closure with respect to $||\cdot||$ coincides to the ideal of compact operators $B_\infty(\cal H)$).
$S(\cal H)$ is closed with respect to $||\cdot ||_1$ and, more strongly, it is a complete metric space with respect to the distance $d_1(\rho,\rho'):= ||\rho-\rho'||_1$.
When $dim(\cal H)$ is finite the three topologies coincide (though the norms do not), as a general result on finite dimensional Banach spaces.
Concerning your last question, there are many viewpoints. My opinion is that a density matrix is physical exactly as pure states are. It is disputable whether or not a mixed state encompasses a sort of physical ignorance, since there is no way to distinguish between "classical probability" and "quantum probability" in a quantum mixture as soon as the mixture is created. See my question Classical and quantum probabilities in density matrices and, in particular Luboš Motl's answer.
See also my answer to Why is the application of probability in QM fundamentally different from application of probability in other areas?
ADDENDUM. In finite dimension, barring the trivial case $dim({\cal H})=2$ where the structure of the space of the states is pictured by the Poincaré-Bloch ball as a manifold with boundary, $S(\cal H)$ has a structure which generalizes that of a manifold with boundary. A stratified space. Roughly speaking, it is not a manifold but is the union of (Riemannian) manifolds with different dimension (depending on the range of the operators) and the intersections are not smooth. When the dimension of $\cal H$ is infinite, one should deal with the notion of infinite dimensional manifold and things become much more complicated.
Let me rephrase your question.
Suppose that $\rho$ is a pure state, i.e., it is written as $$\rho =|\psi\rangle \langle \psi|$$ for some unit vector $\psi$.
Your question is the following.
Q1. Is it possible to find a set of vectors $\phi_1, \ldots, \phi_n$ satisfying $$n>1\:,$$ $||\phi_i||=1$, possibly $\langle \phi_i|\phi_j \rangle \neq 0$ for some $i\neq j$, and numbers $q_1, \ldots, q_n$ with $0< q_i <1$ and $\sum_i q_i =1$, such that
$$|\psi\rangle \langle \psi| = \sum_{i=1}^n q_i|\phi_i \rangle \langle \phi_i|$$
and $|\phi_i \rangle \langle \phi_i| \neq |\phi_j\rangle \langle \phi_j|$
for some $i\neq j$?
As far as I understand you already know the following general result.
THEOREM1. Consider an operator $\rho: H \to H$ where $H$ is a complex Hilbert space and $\rho$ is trace class, non-negative and $tr(\rho)=1$.
Under these hypotheses, $\rho$ is a pure state if and only if $tr(\rho) = tr(\rho^2)$.
As a consequence, since the operator $\sum_{i=1}^n q_i|\phi_i \rangle \langle \phi_i|$ is trace class, non-negative with unit trace, Q1 may be restated as follows.
Q2. Is it possible to find a set of vectors $\phi_1, \ldots, \phi_n$ with $$n>1\:,$$ $||\phi_i||=1$, possibly $\langle \phi_i|\phi_j \rangle \neq 0$ for some $i\neq j$, and numbers $q_1, \ldots, q_n$ with $0< q_i <1$ and $\sum_i q_i =1$, such that
$$tr\left[\left(\sum_{i=1}^n q_i|\phi_i \rangle \langle \phi_i|\right)^2\right] =1$$
and $|\phi_i \rangle \langle \phi_i| \neq |\phi_j\rangle \langle \phi_j|$
for some $i\neq j$?
The answer to Q2 is always negative as soon as $n>1$, and thus
it is not necessary to compute the trace of $\left(\sum_{i=1}^n q_i|\phi_i \rangle \langle \phi_i|\right)^2$, just knowing that $n>1$ is enough to decide that the state $\sum_{i=1}^n q_i|\phi_i \rangle \langle \phi_i|$ is not pure unless $|\phi_i \rangle \langle \phi_i|=|\phi_j \rangle \langle \phi_j|$ for all $i,j$.
The proof is the following. First of all let me introduce the Hilbert-Schmidt scalar product between Hilbert Schmidt operators, and thus trace class operators in particular,
$$(\rho|\rho')_{HS} := tr(\rho^*\rho')\:.$$
The associated norm reads
$$||\rho||_{HS}= \sqrt{tr(\rho^*\rho)}\:.$$
Theorem1 can equivalently be restated as follows.
THEOREM2. Consider an operator $\rho: H \to H$ where $H$ is a complex Hilbert space and $\rho$ is trace class, non-negative and $tr(\rho)=1$.
Under these hypotheses, $\rho$ is a pure state if and only if $||\rho||_{HS}=1$.
Now consider an operator $\rho: H \to H$ of the form
$$\rho = \sum_{i=1}^n q_i \rho_i\tag{0}$$
where $\rho_i := |\phi_i \rangle \langle \phi_i|$ with $\phi_i$ and $q_i$ as in Q2. $\rho$ is trace class, non-negative and we want to check if $||\rho||_{HS}=1$ is possible when $n>1$. This condition is equivalent to saying that $\rho$ is pure.
We can always restrict ourselves to deal with a real vector space of trace class operators, since our trace class operators are self-adjoint and the linear combinations we consider are constructed with real (and non-negative) numbers. The scalar product $(\:|\:)_{HS}$ becomes a standard real (symmetric) scalar product in that real subspace.
The crucial observation is that, as it happens in every real vector space equipped with a real scalar product,
$$\left|\left|\sum_{i=1}^n x_i\right|\right| \leq \sum_{i=1}^n ||x_i||\tag{1}$$
and
"$\leq$" is replaced for "$=$" if and only if $x_i = \alpha_i x$ for some fixed $x$ and non negative numbers $\alpha_i$ where $i=1,\ldots,n$.
In other words,
$$\left|\left|\sum_{i=1}^n q_i \rho_i\right|\right|_{HS} \leq \sum_{i=1}^n ||q_i \rho_i||_{HS}\tag{2}$$
and
"$\leq$" is replaced for "$=$" if and only if $q_i \rho_i = \alpha_i T$ for some fixed $T$ and non negative numbers $\alpha_i$ where $i=1,\ldots,n$.
Since we know that
$$\sum_{i=1}^n ||q_i \rho_i||_{HS}= \sum_{i=1}^n q_i ||\rho_i||_{HS} = \sum_{i=1}^n q_i 1 = \sum_{i=1}^n q_i =1$$
we conclude that If $\rho$ in (0) is pure, then the sign "$\leq$" in (2) is replaced by "$=$", so that $q_i\rho_i = \alpha_i T$ for some fixed operator $T$ and reals $\alpha_i$. Taking the trace of both sides $q_i = \alpha_i tr(T)$ where $tr(T) \neq 0$ because $q_i \neq 0$. Re-defining $T \to \rho_0 := \frac{1}{tr T}T$, we have found that there is a positive trace-class operator $\rho_0$ with unit trace such that $\rho_i= \rho_0$ and furthermore $tr \rho_0^2 = tr \rho_i^2 =1$ so that $\rho_0$ is pure and thus it can be written as $\rho_0 := |\phi_0 \rangle \langle \phi_0|$ for some unit vector $\phi_0$.
Summing up, we have obtained that
if $\rho$ in (0) is pure, then
$|\phi_i\rangle \langle \phi_i|= |\phi_0 \rangle \langle \phi_0|$ for all $i=1,\ldots, n$.
Best Answer
I am not an expert regarding quantum channels etc., but I think this is a pure math question. Let's assume that we're working on a finite-dimensional complex Hilbert space.
In finite dimensions, the defining properties of a density matrix $\rho$ are positive semi-definiteness $\rho \geq 0$ and the trace-normalization $\mathrm{Tr}\,\rho=1$. Properties $1$ and $4$ are consequences of this definition. In particular, a positive semi-definite operator is hermitian$^\dagger$.
To answer your first question: Let $\phi$ denote a map which maps positive semi-definite operators to positive semi-definite operators. Following the considerations above, we see that if $\rho$ is positive semi-definite, then $\phi(\rho)$ is hermitian.
I don't know if I correctly understand your second question, but note that Kraus operators preserve the normalization and positive semi-definiteness. To see this, define for a density operator $\rho$ the following map: $$\varphi(\rho):= \sum\limits_k B_k \,\rho \, B_k^\dagger \quad , $$
with $\displaystyle \sum\limits_k B_k^\dagger\, B_k = \mathbb I$, where $\mathbb I$ denotes the identity operator. We find $$ \mathrm{Tr}\,\varphi(\rho) = \mathrm{Tr}\, \sum\limits_k B_k\, \rho\,B_k^\dagger = \mathrm{Tr}\, \sum\limits_k B_k^\dagger B_k\, \rho = \mathrm{Tr}\, \mathbb I\,\rho = \mathrm{Tr}\,\rho = 1 \quad , $$
where we've used the linearity and the cyclic property of the trace. Moreover:
$$\langle \psi|\varphi(\rho)|\psi\rangle = \sum\limits_k \langle \psi|B_k\, \rho \,B_k^\dagger|\psi\rangle = \sum\limits_k \langle \underbrace{\psi_k|\rho|\psi_k\rangle}_{\geq 0} \geq 0 \quad , $$
with $|\psi_k\rangle := B_k^\dagger|\psi\rangle$.
So $\varphi(\rho)$ is a density matrix again and therefore pure if and only if it is a projector. However, note that Kraus operators can change the purity of a density operator, cf. this quantum computing question.
$^\dagger$ This is easy to see: For a positive semi-definite $\rho$, define the following hermitian operators: \begin{align} \rho_1 &:= \frac{1}{2} (\rho+ \rho^\dagger) \\ \rho_2 &:= \frac{1}{2i} (\rho -\rho^\dagger) \quad . \end{align} Then obviously $\rho = \rho_1 + i\,\rho_2$ and since $\langle \psi |\rho|\psi\rangle , \langle \psi |\rho_1|\psi\rangle , \langle \psi |\rho_2|\psi\rangle \in \mathbb R$, it follows that $\langle \psi|\rho_2|\psi\rangle \overset{!}{=}0$ for all $|\psi\rangle$. Thus $\rho_2=0$ and $\rho=\rho_1$, i.e. $\rho$ is hermitian.