Like many ideas in quantum information theory, I think this is best understood using a $2$-party communication scenario. Suppose Alice has a classical random variable, $X$ which can take values $1,2, \cdots, k$ with probabilities $p_{X}(1), p_{X}(2), \cdots, p_{X}(k)$. Alice then encodes this information by encoding the classical index $j$ in the state $\rho^{j}$. One can represent this scenario as a classical ensemble, $\mathcal{E} = \{ p_X(j), \rho^{j} \}_{j=1}^{k}$ (note that the set $\{\rho^j\}$ is, per se, not mutually orthogonal). For convenience, let's explicitly keep the classical index $j$ by representing this as a classical-quantum state (where the classical index $j$ is correlated to the state $\rho^{j}$ that carries its information) $$ \sigma = \sum\limits_{j=1}^{k} p_X(j) | j \rangle_{X} \langle j | \otimes \rho^{j}. $$

Now, Alice sends this state to Bob, whose task is to determine the classical index $j$ by performing some (optimal) measurement on the state. On some thought, it becomes clear that this is equal to the maximum mutual information of this ensemble. Define, $$ I_{\mathrm{acc}}(\mathcal{E})=\max _{\left\{\Lambda_{y}\right\}} I(X ; Y), $$ where $\{ \Lambda_{y} \}$ is a POVM and $Y$ is a random variable corresponding to the outcome of the measurement. This quantity $I_{\mathrm{acc}}(\mathcal{E})$ is called the accessible information of the ensemble $\mathcal{E}$. Now, in general, one has $$ I_{\mathrm{acc}}(\mathcal{E}) \leq \chi(\mathcal{E}) $$ where $$ \chi(\mathcal{E}) \equiv H\left(\rho_{B}\right)-\sum_{x} p_{X}(x) H\left(\rho_{B}^{x}\right) $$ is the Holevo information --- but this is where out classical-quantum state will become useful. Interestingly, for classical-quantum states, the Holevo information is equal to the mutual information. That is, $$ \chi(\mathcal{E})=I(X ; B)_{\sigma}, $$ which, when combined with the following (simple) bound: $$ I(X:Y) \leq \log \left( \mathrm{dim}(\mathcal{H}) \right), $$ gives us the desired result. Note that the $\mathrm{dim}(\mathcal{H})$ is the Hilbert space where the states $\{\rho^j\}$ belong.

To make the final result transparent, it is instructive to ask what kind of states will saturate this upper bound on the mutual information (and in turn the accessible information). This would correspond to the case where the maximum amount of information can be encoded and accessed from this protocol. It is a simple exercise to show that this happens when the set $\{ \rho_{j} \}$ is mutually orthogonal and hence all states $\rho^j$ are distinguishable. Now, if $k=2^n$, say, for example, because the random variable takes values in $n$-bit strings, then, we need, $\mathrm{dim}(\mathcal{H}) = 2^n$, which can be achieved by choosing states from an $n$-qubit space, $\mathcal{H} \cong (\mathbb{C}^{2})^{\otimes n}$. Therefore, if we want to classically encode (and retrieve) $n$-bits, then we need $n$-qubits. Vice versa, $n$-qubits can contain at most $n$-bits of information.

A few remarks:

1. You don't need to store information in $n$-qubits. You can store the information in any $k$-dimensional quantum system (I'm pointing this out because the tensor product structure of the qubit space plays no role in this protocol, it might as well be a single-particle space with $k$-levels).
2. The key constraint comes from the ability to successfully retrieve information, which requires the states to be distinguishable.

More details can be found in Section 11.6 of Mark Wilde's book.

Like many ideas in quantum information theory, I think this is best understood using a $2$-party communication scenario. Suppose Alice has a classical random variable, $X$ which can take values $1,2, \cdots, k$ with probabilities $p_{X}(1), p_{X}(2), \cdots, p_{X}(k)$. Alice then encodes this information by encoding the classical index $j$ in the state $\rho^{j}$. One can represent this scenario as a classical ensemble, $\mathcal{E} = \{ p_X(j), \rho^{j} \}_{j=1}^{k}$ (note that the set $\{\rho^j\}$ is, per se, not mutually orthogonal). For convenience, let's explicitly keep the classical index $j$ by representing this as a classical-quantum state (where the classical index $j$ is correlated to the state $\rho^{j}$ that carries its information) $$ \sigma = \sum\limits_{j=1}^{k} p_X(j) | j \rangle_{X} \langle j | \otimes \rho^{j}. $$

Now, Alice sends this state to Bob, whose task is to determine the classical index $j$ by performing some (optimal) measurement on the state. On some thought, it becomes clear that this is equal to the maximum mutual information of this ensemble. Define, $$ I_{\mathrm{acc}}(\mathcal{E})=\max _{\left\{\Lambda_{y}\right\}} I(X ; Y), $$ where $\{ \Lambda_{y} \}$ is a POVM and $Y$ is a random variable corresponding to the outcome of the measurement. This quantity $I_{\mathrm{acc}}(\mathcal{E})$ is called the accessible information of the ensemble $\mathcal{E}$. Now, in general, one has $$ I_{\mathrm{acc}}(\mathcal{E}) \leq \chi(\mathcal{E}) $$ where $$ \chi(\mathcal{E}) \equiv H\left(\rho_{B}\right)-\sum_{x} p_{X}(x) H\left(\rho_{B}^{x}\right) $$ is the Holevo information --- but this is where out classical-quantum state will become useful. Interestingly, for classical-quantum states, the Holevo information is equal to the mutual information. That is, $$ \chi(\mathcal{E})=I(X ; B)_{\sigma}, $$ which, when combined with the following (simple) bound: $$ I(X;Y) \leq \log \left( \mathrm{dim}(\mathcal{H}) \right), $$ gives us the desired result. Note that the $\mathrm{dim}(\mathcal{H})$ is the Hilbert space where the states $\{\rho^j\}$ belong.

To make the final result transparent, it is instructive to ask what kind of states will saturate this upper bound on the mutual information (and in turn the accessible information). This would correspond to the case where the maximum amount of information can be encoded and accessed from this protocol. It is a simple exercise to show that this happens when the set $\{ \rho_{j} \}$ is mutually orthogonal and hence all states $\rho^j$ are distinguishable. Now, if $k=2^n$, say, for example, because the random variable takes values in $n$-bit strings, then, we need, $\mathrm{dim}(\mathcal{H}) = 2^n$, which can be achieved by choosing states from an $n$-qubit space, $\mathcal{H} \cong (\mathbb{C}^{2})^{\otimes n}$. Therefore, if we want to classically encode (and retrieve) $n$-bits, then we need $n$-qubits. Vice versa, $n$-qubits can contain at most $n$-bits of information.

A few remarks:

1. You don't need to store information in $n$-qubits. You can store the information in any $k$-dimensional quantum system (I'm pointing this out because the tensor product structure of the qubit space plays no role in this protocol, it might as well be a single-particle space with $k$-levels).
2. The key constraint comes from the ability to successfully retrieve information, which requires the states to be distinguishable.

More details can be found in Section 11.6 of Mark Wilde's book.

Like many ideas in quantum information theory, I think this is best understood using a $2$-party communication scenario. Suppose Alice has a classical random variable, $X$ which can take values $1,2, \cdots, k$ with probabilities $p_{X}(1), p_{X}(2), \cdots, p_{X}(k)$. Alice then encodes this information by encoding the classical index $j$ in the state $\rho^{j}$. One can represent this scenario as a classical ensemble, $\mathcal{E} = \{ p_X(j), \rho^{j} \}_{j=1}^{k}$ (note that the set $\{\rho^j\}$ is, per se, not mutually orthogonal). For convenience, let's explicitly keep the classical index $j$ by representing this as a classical-quantum state (where the classical index $j$ is correlated to the state $\rho^{j}$ that carries its information) $$ \sigma = \sum\limits_{j=1}^{k} p_X(j) | j \rangle_{X} \langle j | \otimes \rho^{j}. $$

Now, Alice sends this state to Bob, whose task is to determine the classical index $j$ by performing some (optimal) measurement on the state. On some thought, it becomes clear that this is equal to the maximum mutual information of this ensemble. Define, $$ I_{\mathrm{acc}}(\mathcal{E})=\max _{\left\{\Lambda_{y}\right\}} I(X ; Y), $$ where $\{ \Lambda_{y} \}$ is a POVM and $Y$ is a random variable corresponding to the outcome of the measurement. This quantity $I_{\mathrm{acc}}(\mathcal{E})$ is called the accessible information of the ensemble $\mathcal{E}$. Now, in general, one has $$ I_{\mathrm{acc}}(\mathcal{E}) \leq \chi(\mathcal{E}) $$ where $$ \chi(\mathcal{E}) \equiv H\left(\rho_{B}\right)-\sum_{x} p_{X}(x) H\left(\rho_{B}^{x}\right) $$ is the Holevo information --- but this is where out classical-quantum state will become useful. Interestingly, for classical-quantum states, the Holevo information is equal to the mutual information. That is, $$ \chi(\mathcal{E})=I(X ; B)_{\sigma}, $$ which, when combined with the following (simple) bound: $$ I(X:Y) \leq \log \left( \mathrm{dim}(\mathcal{H}) \right), $$ gives us the desired result.

To make the final result transparent, it is instructive to ask what kind of states will saturate this upper bound on the mutual information (and in turn the accessible information). This would correspond to the case where the maximum amount of information can be encoded and accessed from this protocol. It is a simple exercise to show that this happens when the set $\{ \rho_{j} \}$ is mutually orthogonal and hence all states $\rho^j$ are distinguishable, maximum number of which are $2^{n}$ when we are encoding this information in an $n$-qubit space --- thereby proving the claim that $n$-qubits can contain at most $n$-bits of information.

More details can be found in Section 11.6 of Mark Wilde's book.

Like many ideas in quantum information theory, I think this is best understood using a $2$-party communication scenario. Suppose Alice has a classical random variable, $X$ which can take values $1,2, \cdots, k$ with probabilities $p_{X}(1), p_{X}(2), \cdots, p_{X}(k)$. Alice then encodes this information by encoding the classical index $j$ in the state $\rho^{j}$. One can represent this scenario as a classical ensemble, $\mathcal{E} = \{ p_X(j), \rho^{j} \}_{j=1}^{k}$ (note that the set $\{\rho^j\}$ is, per se, not mutually orthogonal). For convenience, let's explicitly keep the classical index $j$ by representing this as a classical-quantum state (where the classical index $j$ is correlated to the state $\rho^{j}$ that carries its information) $$ \sigma = \sum\limits_{j=1}^{k} p_X(j) | j \rangle_{X} \langle j | \otimes \rho^{j}. $$

Now, Alice sends this state to Bob, whose task is to determine the classical index $j$ by performing some (optimal) measurement on the state. On some thought, it becomes clear that this is equal to the maximum mutual information of this ensemble. Define, $$ I_{\mathrm{acc}}(\mathcal{E})=\max _{\left\{\Lambda_{y}\right\}} I(X ; Y), $$ where $\{ \Lambda_{y} \}$ is a POVM and $Y$ is a random variable corresponding to the outcome of the measurement. This quantity $I_{\mathrm{acc}}(\mathcal{E})$ is called the accessible information of the ensemble $\mathcal{E}$. Now, in general, one has $$ I_{\mathrm{acc}}(\mathcal{E}) \leq \chi(\mathcal{E}) $$ where $$ \chi(\mathcal{E}) \equiv H\left(\rho_{B}\right)-\sum_{x} p_{X}(x) H\left(\rho_{B}^{x}\right) $$ is the Holevo information --- but this is where out classical-quantum state will become useful. Interestingly, for classical-quantum states, the Holevo information is equal to the mutual information. That is, $$ \chi(\mathcal{E})=I(X ; B)_{\sigma}, $$ which, when combined with the following (simple) bound: $$ I(X:Y) \leq \log \left( \mathrm{dim}(\mathcal{H}) \right), $$ gives us the desired result. Note that the $\mathrm{dim}(\mathcal{H})$ is the Hilbert space where the states $\{\rho^j\}$ belong.

To make the final result transparent, it is instructive to ask what kind of states will saturate this upper bound on the mutual information (and in turn the accessible information). This would correspond to the case where the maximum amount of information can be encoded and accessed from this protocol. It is a simple exercise to show that this happens when the set $\{ \rho_{j} \}$ is mutually orthogonal and hence all states $\rho^j$ are distinguishable. Now, if $k=2^n$, say, for example, because the random variable takes values in $n$-bit strings, then, we need, $\mathrm{dim}(\mathcal{H}) = 2^n$, which can be achieved by choosing states from an $n$-qubit space, $\mathcal{H} \cong (\mathbb{C}^{2})^{\otimes n}$. Therefore, if we want to classically encode (and retrieve) $n$-bits, then we need $n$-qubits. Vice versa, $n$-qubits can contain at most $n$-bits of information.

A few remarks:

1. You don't need to store information in $n$-qubits. You can store the information in any $k$-dimensional quantum system (I'm pointing this out because the tensor product structure of the qubit space plays no role in this protocol, it might as well be a single-particle space with $k$-levels).
2. The key constraint comes from the ability to successfully retrieve information, which requires the states to be distinguishable.

More details can be found in Section 11.6 of Mark Wilde's book.