text:
# The capacity of the quantum depolarizing channel

## 1 Background and statement of results

### 1.1 Introduction

This paper computes the capacity of the $`d`$-dimensional quantum depolarizing channel for transmission of classical information. The result confirms a longstanding conjecture, namely that the best rate of information transfer can be achieved without any entanglement across multiple uses of the channel. It is sufficient to choose an orthonormal basis for the state space, and then use this basis both to encode messages as product states at the input side and to perform measurements at the output side which project onto this same basis. In this sense the depolarizing channel can be treated as a classical channel.

The Holevo-Schumacher-Westmoreland theorem allows the capacity to be expressed in terms of the Holevo quantity $`\chi ^{}`$. If the Holevo quantity is additive, then this expression implies that the capacity is actually equal to $`\chi ^{}`$. In this paper we prove several additivity properties for the depolarizing channel, including the additivity of $`\chi ^{}`$. One important mathematical tool used in the proof is the Lieb-Thirring inequality, which provides a bound for the non-commutative $`p`$-norm of a product of positive matrices. These notions, as well as the definition of channel capacity and its relation to the Holevo quantity, are described in the following subsections.

### 1.2 The depolarizing channel

The depolarizing channel is a particularly simple model for noise in quantum systems , and has been studied in a variety of contexts . In $`d`$ dimensions the model is implemented by a completely positive trace-preserving map $`\mathrm{\Delta }_\lambda `$, depending on one real parameter $`\lambda `$, which maps a state $`\rho `$ on $`𝐂^d`$ into a linear combination of itself and the $`d\times d`$ identity matrix $`I`$:

$`\mathrm{\Delta }_\lambda (\rho )=\lambda \rho +{\displaystyle \frac{1\lambda }{d}}I`$ (1)
The condition of complete positivity requires that $`\lambda `$ satisfy the bounds

$`{\displaystyle \frac{1}{d^21}}\lambda 1`$ (2)
The channel $`\mathrm{\Delta }_\lambda `$ maps a pure input state to a mixed output state. Because the channel is highly symmetric, all such output states are unitarily equivalent, and have eigenvalues $`\lambda +(1\lambda )/d`$ (with multiplicity 1) and $`(1\lambda )/d`$ (with multiplicity $`d1`$).

### 1.3 Measures of noisiness

We will use three measures of noisiness for quantum channels, namely the minimal output entropy $`S_{\mathrm{min}}`$, the maximal output $`p`$-norm $`\nu _p`$, and the Holevo quantity $`\chi ^{}`$. First recall that a state $`\rho `$ is a positive operator with trace equal to $`1`$, and its von Neumann entropy is defined by

$`S(\rho )=\mathrm{Tr}\rho \mathrm{log}\rho `$ (3)
The minimal output entropy of the channel $`\mathrm{\Psi }`$ is defined as

$`S_{\mathrm{min}}(\mathrm{\Psi })=\underset{\rho }{inf}S\left(\mathrm{\Psi }(\rho )\right)`$ (4)
Recall also that the $`p`$-norm of a positive matrix $`A`$ is defined for $`p1`$ by

$`A_p=\left(\mathrm{Tr}A^p\right)^{1/p}`$ (5)
Amosov, Holevo and Werner  introduced the notion of the maximal $`p`$-norm of a channel as a way to characterize its noisiness. This quantity is defined as

$`\nu _p(\mathrm{\Psi })=\underset{\rho }{sup}\mathrm{\Psi }(\rho )_p`$ (6)
Since the entropy of a state is the negative of the derivative of the $`p`$-norm at $`p=1`$, it follows that for any channel $`\mathrm{\Psi }`$,

$`{\displaystyle \frac{d}{dp}}\nu _p(\mathrm{\Psi })|_{p=1}=S_{\mathrm{min}}(\mathrm{\Psi })`$ (7)
The third measure of noisiness, the Holevo quantity, is closely related to the information-carrying capacity of the channel, as will be explained in the next section. We will use the symbol $``$ to denote an ensemble of input states for the channel, that is a collection of states $`\rho _i`$ together with a probability distribution $`\pi _i`$. The Holevo quantity is

$`\chi ^{}(\mathrm{\Psi })=\underset{}{sup}\left[S\left(\mathrm{\Psi }(\rho )\right){\displaystyle \underset{i}{}}\pi _iS\left(\mathrm{\Psi }(\rho _i)\right)\right]`$ (8)
where $`\rho =\pi _i\rho _i`$ is the average input state of the ensemble.

These three measures can be easily computed for the depolarizing channel. The values are

$`S_{\mathrm{min}}(\mathrm{\Delta }_\lambda )=\left(\lambda +{\displaystyle \frac{1\lambda }{d}}\right)\mathrm{log}\left(\lambda +{\displaystyle \frac{1\lambda }{d}}\right)(d1)({\displaystyle \frac{1\lambda }{d}})\mathrm{log}\left({\displaystyle \frac{1\lambda }{d}}\right)`$ (9)
$`\nu _p(\mathrm{\Delta }_\lambda )=\left[\left(\lambda +{\displaystyle \frac{1\lambda }{d}}\right)^p+(d1)\left({\displaystyle \frac{1\lambda }{d}}\right)^p\right]^{1/p}`$ (10)
$`\chi ^{}(\mathrm{\Delta }_\lambda )=\mathrm{log}dS_{\mathrm{min}}(\mathrm{\Delta }_\lambda )`$ (11)
The value (11) is achieved by choosing an ensemble $``$ consisting of pure states belonging to an orthonormal basis (it does’nt matter which one), and choosing the uniform distribution $`\pi _i=1/d`$. Since the average input state for this ensemble is $`(1/d)I`$, the terms on the right side of (8) are separately maximized for this choice of ensemble, and this leads to the result (11).

### 1.4 Additivity conjectures

It is conjectured that $`S_{\mathrm{min}}`$ and $`\chi ^{}`$ are both additive for product channels. This would mean that for any channels $`\mathrm{\Psi }_1`$ and $`\mathrm{\Psi }_2`$,

$`S_{\mathrm{min}}(\mathrm{\Psi }_1\mathrm{\Psi }_2)=S_{\mathrm{min}}(\mathrm{\Psi }_1)+S_{\mathrm{min}}(\mathrm{\Psi }_2)`$ (12)
and

$`\chi ^{}(\mathrm{\Psi }_1\mathrm{\Psi }_2)=\chi ^{}(\mathrm{\Psi }_1)+\chi ^{}(\mathrm{\Psi }_2)`$ (13)
Equivalently, the conjectures would imply that for product channels both of these measures of noisiness are achieved with product input states. As we explain in the next section, the Holevo quantity $`\chi ^{}`$ is related to channel capacity, and in fact (13) would imply that the channel capacity of an arbitrary channel $`\mathrm{\Psi }`$ is precisely $`\chi ^{}(\mathrm{\Psi })`$. These conjectures have been established for some special classes of channels, including all unital qubit channels  and all entanglement-breaking channels . However a general proof has remained elusive.

It was further conjectured in  that the quantity $`\nu _p`$ should be multiplicative for product channels for all $`p1`$, reflecting the idea that the $`p`$-norm of a product channel would be maximized with product states. This conjecture would imply (12), since the entropy can be obtained from the derivative of the $`p`$-norm at $`p=1`$. Again this conjecture has been established for some special classes of channels, in particular for unital qubit channels . However it is now known that the conjecture does not hold in general – this was demonstrated recently by the discovery of a family of counterexamples for values $`p5`$ . Nevertheless, in this paper we will prove that the AHW conjecture is true for any product channel of the form $`\mathrm{\Delta }_\lambda \mathrm{\Psi }`$ where $`\mathrm{\Psi }`$ is arbitrary, and we will show how this property implies the additivity of $`S_{\mathrm{min}}`$ and of $`\chi ^{}`$ for such a product channel.

### 1.5 Channel capacity

In order to relate these additivity results to the channel capacity problem, we recall first the definition of the capacity for a general channel $`\mathrm{\Psi }`$. Again we denote by $``$ an ensemble of input states, and we also denote by $``$ a measurement, or POVM, at the output side of the channel. Recall that a POVM is a collection of positive operators $`\{E_j\}`$ which sum to the identity matrix. When a state $`\rho `$ is measured using a POVM $`\{E_j\}`$, the outcome $`j`$ is obtained with probability $`\mathrm{Tr}(\rho E_j)`$. This notion generalizes the familiar von Neumann measurement, which is the special case where the operators $`E_j`$ are orthogonal projections.

The ensemble $``$, the POVM $``$ and the channel $`\mathrm{\Psi }`$ together define a classical noisy channel, whose transition matrix is

$`p_{ij}=\mathrm{Tr}\left[\mathrm{\Psi }(\rho _i)E_j\right]`$ (14)
If we write $`X`$ for a random input signal with distribution $`\pi _i=P(X=i)`$, then the output signal $`Y`$ from this classical channel has distribution

$`P(Y=j)={\displaystyle \underset{i}{}}\pi _ip_{ij}`$ (15)
The Shannon capacity of a classical noisy channel measures the maximum rate at which information can be reliably transmitted through the channel. Shannon’s formula computes this capacity as the maximum of the mutual information $`I(X,Y)`$ between an input signal $`X`$ and its corresponding output signal $`Y`$ given by (15), where the maximum is evaluated over all choices of distribution $`\{\pi _i\}`$ for $`X`$. For the case of a quantum channel, we are interested in the maximum rate that can be achieved using the optimal choices of input states $`\{\rho _i\}`$ and of output measurements $`\{E_j\}`$. Therefore we are led to define the Shannon capacity of the quantum channel $`\mathrm{\Psi }`$ as

$`C_{\mathrm{Shan}}(\mathrm{\Psi })=\underset{,}{sup}I(X,Y)`$ (16)
where the distribution of the input signal $`X`$ is determined by the ensemble $``$, that is $`P(X=i)=\pi _i`$, and where the output signal $`Y`$ is determined by (15). This can be easily evaluated for the depolarizing channel. The mutual imformation $`I(X,Y)`$ in (16) is maximized by choosing an ensemble consisting of projections onto an orthonormal basis, and using the same basis for the measurement. The result is

$`C_{\mathrm{Shan}}(\mathrm{\Delta }_\lambda )=\chi ^{}(\mathrm{\Delta }_\lambda )`$ (17)
where the capacity $`\chi ^{}(\mathrm{\Delta }_\lambda )`$ is given by (11).

If two copies of the channel $`\mathrm{\Psi }`$ are available then it may be possible to achieve a higher rate of transmission by sharing the signals across the two channels. This possibility exists because quantum channels have an additional resource which is not available for classical channels, namely entangled states which can be used to encode signals for the product channel. It is also possible to make measurements at the output side using a POVM which projects onto entangled states. With these resources the best rate that can be achieved using two copies of the channel is

$`{\displaystyle \frac{1}{2}}C_{\mathrm{Shan}}(\mathrm{\Psi }\mathrm{\Psi })`$ (18)
It is known that in general (18) is larger than (16) . This observation leads to the question of finding the asympototic capacity which would be achieved by sharing the input signals across an unlimited number of copies of $`\mathrm{\Psi }`$. This ultimate capacity is given by

$`C_{\mathrm{ult}}(\mathrm{\Psi })=\underset{n\mathrm{}}{lim}{\displaystyle \frac{1}{n}}C_{\mathrm{Shan}}(\mathrm{\Psi }^n)`$ (19)
(a standard subadditivity argument shows the existence of this limit).

At the present time it is an open problem to determine $`C_{\mathrm{ult}}(\mathrm{\Psi })`$ for an arbitrary channel $`\mathrm{\Psi }`$. However it can be expressed in terms of the Holevo quantity (8). Recall that the Holevo bound implies that

$`C_{\mathrm{Shan}}(\mathrm{\Psi })\chi ^{}(\mathrm{\Psi }),`$ (20)
and hence that

$`C_{\mathrm{ult}}(\mathrm{\Psi })\underset{n\mathrm{}}{lim}{\displaystyle \frac{1}{n}}\chi ^{}(\mathrm{\Psi }^n)`$ (21)
Furthermore the Holevo-Schumacher-Westmoreland theorem  shows that the rate $`\chi ^{}(\mathrm{\Psi })`$ can be achieved with multiple copies of the channel, by restricting to product states for the input signals, but allowing entangled measurements at the outputs. Applying this theorem to the product channel $`\mathrm{\Psi }^n`$ implies that the rate $`(1/n)\chi ^{}(\mathrm{\Psi }^n)`$ is achieved with input signals which may be entangled across $`n`$ uses of the channel. Allowing $`n`$ to be arbitrarily large leads to the result

$`C_{\mathrm{ult}}(\mathrm{\Psi })=\underset{n\mathrm{}}{lim}{\displaystyle \frac{1}{n}}\chi ^{}(\mathrm{\Psi }^n)`$ (22)
### 1.6 Statement of results

Our first result is the evaluation of (19) for the depolarizing channel.

###### Theorem 1

The capacity of the $`d`$-dimensional depolarizing channel is

$`C_{\mathrm{ult}}(\mathrm{\Delta }_\lambda )=\chi ^{}(\mathrm{\Delta }_\lambda )=C_{\mathrm{Shan}}(\mathrm{\Delta }_\lambda )=\mathrm{log}dS_{\mathrm{min}}(\mathrm{\Delta }_\lambda )`$ (23)
where $`S_{\mathrm{min}}(\mathrm{\Delta }_\lambda )`$ is evaluated in (9).

The fact that $`C_{\mathrm{ult}}(\mathrm{\Delta }_\lambda )=C_{\mathrm{Shan}}(\mathrm{\Delta }_\lambda )`$ means that as far as the information-carrying properties of the depolarizing channel are concerned, there is no advantage gained by using either entangled input states or using entangled measurements. The optimal rate can be achieved by choosing an orthonormal basis (because $`\mathrm{\Delta }_\lambda `$ is symmetric it does’nt matter which one) and using this basis to encode the signals and also to measure them. In this sense the channel behaves like a classical channel, and entanglement does not play any role in its capacity.

The basic ingredient in the proof of Theorem 1 is the additivity of the Holevo quantity $`\chi ^{}`$ for the depolarizing channel, which we state in the next Theorem.

###### Theorem 2

For any channel $`\mathrm{\Psi }`$,

$`\chi ^{}(\mathrm{\Delta }_\lambda \mathrm{\Psi })=\chi ^{}(\mathrm{\Delta }_\lambda )+\chi ^{}(\mathrm{\Psi })`$ (24)
Theorem 1 follows easily from this, as we now demonstrate. The result $`C_{\mathrm{ult}}(\mathrm{\Delta }_\lambda )=\chi ^{}(\mathrm{\Delta }_\lambda )`$ in Theorem 1 follows immediately from Theorem 2 by choosing $`\mathrm{\Psi }=\mathrm{\Delta }_\lambda ^n`$ in (24) and applying (22). The second equality $`\chi ^{}(\mathrm{\Delta }_\lambda )=C_{\mathrm{Shan}}(\mathrm{\Delta }_\lambda )`$ was derived in (17).

Finally we state the AHW conjecture for the depolarizing channel, which underlies all the other results. Since the derivative of $`\nu _p(\mathrm{\Psi })`$ at $`p=1`$ is equal to $`S_{\mathrm{min}}(\mathrm{\Psi })`$, the additivity of $`S_{\mathrm{min}}`$ is a special case of the AHW conjecture. We state both results next in Theorem 3.

###### Theorem 3

For any channel $`\mathrm{\Psi }`$, and any $`p1`$,

$`\nu _p(\mathrm{\Delta }_\lambda \mathrm{\Psi })=\nu _p(\mathrm{\Delta }_\lambda )\nu _p(\mathrm{\Psi })`$ (25)
and hence

$`S_{\mathrm{min}}(\mathrm{\Delta }_\lambda \mathrm{\Psi })=S_{\mathrm{min}}(\mathrm{\Delta }_\lambda )+S_{\mathrm{min}}(\mathrm{\Psi })`$ (26)
Special cases of Theorem 3 were previously established, namely for integer values of $`p`$ in all dimensions $`d`$ , and for all $`p1`$ in dimension $`d=2`$ .

### 1.7 Organization

The paper is organised as follows. Section 2 outlines the proof of Theorem 3, and states two key results which are used, namely the convex decomposition of the depolarizing channel, and the bound for the phase-damping channel. These results are then established in Sections 3 and 4, and finally Theorem 2 is proved in Section 5. Section 3 also describes in detail the convex decompositions for the two-dimensional qubit depolarizing channel. Section 6 contains some discussion of the nature of the proof, and some ideas about further directions to pursue.

## 2 Outline of the proof

### 2.1 Definition of phase-damping channel

As discussed above, Theorem 1 follows immediately from Theorem 2, using the HSW Theorem (22). Theorem 2 itself is a slight extension of Theorem 3, and will be proved in Section 5. Most of the work in this paper goes into the proof of the AHW conjecture (25) in Theorem 3. The proof presented here develops further the methods introduced in  where the same result was established for unital qubit channels. The basic idea is similar: we express the depolarizing channel $`\mathrm{\Delta }_\lambda `$ as a convex combination of simpler channels, and then we prove a bound for these simpler channels which implies the result (25). In  the simpler channels were unitarily equivalent to phase-damping channels, and we use the same name for the channels here, which are defined as follows.

###### Definition 4

Let $`=\{|\psi _i\}`$ be an orthonormal basis, and let $`E_i=|\psi _i\psi _i|`$. The phase-damping channel corresponding to $``$ is the one-parameter family of maps

$`\mathrm{\Phi }_\lambda (\rho )=\lambda \rho +(1\lambda ){\displaystyle \underset{i=1}{\overset{d}{}}}E_i\rho E_i`$ (27)
where the parameter $`\lambda `$ satisfies the bounds

$`{\displaystyle \frac{1}{d1}}\lambda 1`$ (28)
The parameter range (28) is required by the condition of complete positivity. If we write $`\rho =(\rho _{ij})`$ as a matrix in the basis $`|\psi _i`$ then the channel (27) acts by scaling the off-diagonal entries and leaving unchanged the diagonal entries, that is

$`\mathrm{\Phi }_\lambda (\rho )_{ij}=\{\begin{array}{cc}\rho _{ij}\hfill & \text{if }i=j\hfill \\ \lambda \rho _{ij}\hfill & \text{if }ij\hfill \end{array}`$ (29)
We will express $`\mathrm{\Delta }_\lambda `$ as a convex combination of phase-damping channels, all with the same parameter $`\lambda `$. Furthermore these phase-damping channels will all share a common property, which is expressed by the following definition.

###### Definition 5

We say that a vector $`𝐯=(v_1,\mathrm{},v_d)`$ in $`𝐂^d`$ is uniform if $`|v_i|=|v_j|`$ for all $`i,j=1,\mathrm{},d`$.

It will turn out that the phase-damping channels which arise in the convex decomposition are constructed from orthonormal bases $`=\{|\psi _i\}`$ where all the vectors $`|\psi _i`$ are uniform. Since each vector $`|\psi _i`$ is normalized, it follows that all its entries have absolute value $`1/\sqrt{d}`$. As a consequence, if $`D`$ is any diagonal matrix, then for any uniform state $`|\psi `$

$`\mathrm{Tr}\left[|\psi \psi |D\right]=\psi |D|\psi ={\displaystyle \frac{1}{d}}\mathrm{Tr}D`$ (30)
###### Definition 6

Let $`\mathrm{\Phi }`$ be the phase-damping channel corresponding to the orthonormal basis $``$. We say that $`\mathrm{\Phi }`$ is uniform if $`|\psi _i`$ is uniform for every $`|\psi _i`$.

### 2.2 Three lemmas

There are three steps in the proof of Theorem 3. The goal is to find a bound for $`(\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})_p`$ which will lead to (25), where $`\mathrm{\Psi }`$ is any other channel, and $`\tau _{12}`$ is any state. The first step is a partial diagonalization of the state $`\tau _{12}`$. This step uses the following invariance property of the depolarizing channel.

###### Lemma 7

Let $`\tau _1=\mathrm{Tr}_2(\tau _{12})`$ denote the reduced density matrix of $`\tau _{12}`$, and let $`U`$ be a unitary matrix. Define $`\tau _{12}^{}=(UI)\tau _{12}(U^{}I)`$. Then for all $`p1`$

$`(\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})_p=(\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12}^{})_p`$ (31)
Proof: the definition of $`\mathrm{\Delta }_\lambda `$ in (1) implies that the unitary matrix $`UI`$ can be pulled through the channel $`\mathrm{\Delta }_\lambda \mathrm{\Psi }`$, and then the invariance of the $`p`$-norm implies (31). QED

The second step uses the following result which expresses $`\mathrm{\Delta }_\lambda `$ as a convex combination of phase-damping channels. This result will be derived in Section 3.

###### Lemma 8

For $`n=1,\mathrm{},2d^2(d+1)`$, there are positive numbers $`c_n`$, unitary matrices $`U_n`$ and uniform phase-damping channels $`\mathrm{\Phi }_\lambda ^{(n)}`$ such that for any state $`\rho `$

$`\mathrm{\Delta }_\lambda (\rho )={\displaystyle \underset{n=1}{\overset{2d^2(d+1)}{}}}c_nU_n^{}\mathrm{\Phi }_{\lambda }^{}{}_{}{}^{(n)}(\rho )U_n`$ (32)
The third step in the proof uses the following bound for the phase-damping channels. This bound will be derived in Section 4.

###### Lemma 9

Let $`\mathrm{\Phi }_\lambda `$ be a phase-damping channel defined as in (27), with corresponding orthogonal projectors $`E_i=|\psi _i\psi _i|`$. For an arbitrary bi-partite state $`\rho _{12}`$ define

$`\rho _2^{(i)}=\mathrm{Tr}_1\left[(E_iI)\rho _{12}\right]`$ (33)
where $`\mathrm{Tr}_1`$ is the trace over the first factor. Recall the factor $`\nu _p(\mathrm{\Delta }_\lambda )`$ from (10). Then for all $`p1`$,

$`(\mathrm{\Phi }_\lambda I)(\rho _{12})_pd^{(11/p)}\nu _p(\mathrm{\Delta }_\lambda )\left[{\displaystyle \underset{i=1}{\overset{d}{}}}\mathrm{Tr}\left(\rho _2^{(i)}\right)^p\right]^{1/p}`$ (34)
### 2.3 Proof of Theorem 3

We will now prove Theorem 3 using these three lemmas. Since the left side of (25) is at least as big as the right side, it is sufficient to prove that for any bipartite state $`\tau _{12}`$

$`(\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})_p\nu _p(\mathrm{\Delta }_\lambda )\nu _p(\mathrm{\Psi })`$ (35)
The first step is to use Lemma 1 to partially diagonalize the state $`\tau _{12}`$. Let $`U`$ be a unitary matrix which diagonalizes the reduced density matrix $`\tau _1=\mathrm{Tr}_2(\tau _{12})`$, and let $`\tau _{12}^{}=(UI)\tau _{12}(U^{}I)`$, so that $`\tau _1^{}=U\tau _1U^{}`$ is diagonal. By Lemma 7 we can replace $`\tau _{12}`$ by $`\tau _{12}^{}`$ without changing the left side of (35). Therefore we will assume henceforth without loss of generality that $`\tau _1`$ is diagonal.

The second step is to apply the convex decomposition (32) on the left side of (35):

$`(\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})={\displaystyle \underset{n=1}{\overset{2d^2(d+1)}{}}}c_n(U_n^{}I)(\mathrm{\Phi }_{\lambda }^{}{}_{}{}^{(n)}\mathrm{\Psi })(\tau _{12})(U_nI)`$ (36)
For the third step, notice that by convexity of the $`p`$-norm it is sufficient to prove the bound (35) for each term $`(\mathrm{\Phi }_{\lambda }^{}{}_{}{}^{(n)}\mathrm{\Psi })(\tau _{12})`$ appearing on the right side of (36), namely

$`(\mathrm{\Phi }_\lambda ^{(n)}\mathrm{\Psi })(\tau _{12})_p\nu _p(\mathrm{\Delta }_\lambda )\nu _p(\mathrm{\Psi })`$ (37)
In order to derive (37), we apply (34) with

$`\rho _{12}=(I\mathrm{\Psi })(\tau _{12}),\rho _2^{(i)}=\mathrm{\Psi }(\tau _2^{(i)})=\mathrm{\Psi }\left(\mathrm{Tr}_1\left[(E_iI)\tau _{12}\right]\right)`$ (38)
Therefore (34) gives

$`(\mathrm{\Phi }_\lambda ^{(n)}\mathrm{\Psi })(\tau _{12})_pd^{(11/p)}\nu _p(\mathrm{\Delta }_\lambda )\left[{\displaystyle \underset{i=1}{\overset{d}{}}}\mathrm{Tr}\left(\mathrm{\Psi }(\tau _2^{(i)})\right)^p\right]^{1/p}`$ (39)
Now the definition of the $`p`$-norm $`\nu _p(\mathrm{\Psi })`$ implies that for each $`i`$,

$`\left[\mathrm{Tr}\left(\mathrm{\Psi }(\tau _2^{(i)})\right)^p\right]^{1/p}\nu _p(\mathrm{\Psi })\mathrm{Tr}(\tau _2^{(i)})`$ (40)
From the definition of $`\tau _2^{(i)}`$ it follows that

$`\mathrm{Tr}(\tau _2^{(i)})=\mathrm{Tr}(E_i\tau _1)`$ (41)
Furthermore in the first step we chose the state $`\tau _1`$ to be diagonal, and from Lemma 8 the phase-damping channels appearing on the right side of (36) are all uniform (recall Definition 6). Hence from (30) we get

$`\mathrm{Tr}(E_i\tau _1)={\displaystyle \frac{1}{d}}\mathrm{Tr}\tau _1={\displaystyle \frac{1}{d}}`$ (42)
Inserting into (39) gives

$`(\mathrm{\Phi }_\lambda ^{(n)}\mathrm{\Psi })(\tau _{12})_p`$ $``$ $`d^{(11/p)}\nu _p(\mathrm{\Delta }_\lambda )\nu _p(\mathrm{\Psi })\left[d\left({\displaystyle \frac{1}{d}}\right)^p\right]^{1/p}`$ (43)
$`=`$ $`\nu _p(\mathrm{\Delta }_\lambda )\nu _p(\mathrm{\Psi })`$ (44)
which completes the proof. QED

## 3 The convex decomposition

### 3.1 Proof of Lemma 9

The derivation proceeds in two stages. First we define a new channel $`\mathrm{\Omega }_\lambda `$ which appears in an intermediate role:

$`\mathrm{\Omega }_\lambda (\rho )=\mathrm{\Delta }_\lambda (\rho )+{\displaystyle \frac{1\lambda }{d}}\left[\rho \mathrm{diag}(\rho )\right]`$ (45)
where $`\mathrm{diag}(\rho )`$ is the diagonal part of the matrix $`\rho `$. It is not hard to see that $`\mathrm{\Omega }_\lambda `$ is completely positive and trace-preserving for all $`\lambda `$ in the range (28). For example it can be re-written as

$`\mathrm{\Omega }_\lambda (\rho )=\left(\lambda +{\displaystyle \frac{1\lambda }{d}}\right)\rho +{\displaystyle \frac{(d1)(1\lambda )}{d}}{\displaystyle \frac{1}{d1}}\left[I\mathrm{diag}(\rho )\right],`$ (46)
and the map $`\left[I\mathrm{diag}(\rho )\right]/(d1)`$ is easily seen to be completely positive and trace-preserving. Next let $`G`$ be the diagonal unitary matrix with entries

$`G_{kk}`$ $`=`$ $`\mathrm{exp}\left({\displaystyle \frac{2\pi ik}{d}}\right),1kd`$ (47)
$`G_{kl}`$ $`=`$ $`0,kl`$ (48)
###### Lemma 10

For any matrix $`\rho `$,

$`\mathrm{\Delta }_\lambda (\rho )={\displaystyle \frac{\lambda d}{1+(d1)\lambda }}\mathrm{\Omega }_\lambda (\rho )+{\displaystyle \frac{1\lambda }{1+(d1)\lambda }}{\displaystyle \frac{1}{d}}{\displaystyle \underset{k=1}{\overset{d}{}}}\left(G^{}\right)^k\mathrm{\Omega }_\lambda (\rho )G^k`$ (49)
Proof: A straightforward computation shows that for any matrix $`\rho `$

$`{\displaystyle \frac{1}{d}}{\displaystyle \underset{k=1}{\overset{d}{}}}\left(G^{}\right)^k\rho G^k=\mathrm{diag}(\rho )`$ (50)
Applying the definitions of $`\mathrm{\Delta }_\lambda `$ and $`\mathrm{\Omega }_\lambda `$ from (1) and (45), the result now follows easily. QED

For the second stage in the derivation of (32) we express the channel $`\mathrm{\Omega }_\lambda `$ itself as a convex combination of phase-damping channels. To this end define the diagonal unitary matrix $`H`$ by

$`H_{kk}`$ $`=`$ $`\mathrm{exp}\left({\displaystyle \frac{2\pi ik^2}{2d^2}}\right),1kd`$ (51)
$`H_{kl}`$ $`=`$ $`0,kl`$ (52)
and define the following pure state $`|\theta `$:

$`|\theta ={\displaystyle \frac{1}{\sqrt{d}}}\left(\begin{array}{c}1\\ 1\\ \mathrm{}\\ 1\end{array}\right)`$ (53)
For each $`k=1,\mathrm{},d`$ and $`a=1,\mathrm{},2d^2`$ we define the pure state

$`|\psi _{k,a}=G^kH^a|\theta `$ (54)
and the corresponding orthogonal projection

$`E_{k,a}=|\psi _{k,a}\psi _{k,a}|`$ (55)
For each fixed $`a`$, the states $`\{|\psi _{k,a}\}`$ form an orthonormal basis. We denote by $`\mathrm{\Phi }_\lambda ^{(a)}`$ the corresponding family of phase-damping channels, that is

$`\mathrm{\Phi }_\lambda ^{(a)}(\rho )=\lambda \rho +(1\lambda ){\displaystyle \underset{k=1}{\overset{d}{}}}E_{k,a}\rho E_{k,a},a=1,\mathrm{},2d^2`$ (56)
###### Lemma 11

$`\mathrm{\Omega }_\lambda ={\displaystyle \frac{1}{2d^2}}{\displaystyle \underset{a=1}{\overset{2d^2}{}}}\mathrm{\Phi }_\lambda ^{(a)}`$ (57)
Proof: Using the definitions of $`\mathrm{\Omega }_\lambda `$, $`\mathrm{\Delta }_\lambda `$ and $`\mathrm{\Phi }_\lambda ^{(a)}`$, it suffices to show that for any state $`\rho `$

$`{\displaystyle \frac{1}{2d}}{\displaystyle \underset{a=1}{\overset{2d^2}{}}}{\displaystyle \underset{k=1}{\overset{d}{}}}E_{k,a}\rho E_{k,a}=I+\rho \mathrm{diag}(\rho )`$ (58)
For each $`x=1,\mathrm{},d`$, let $`|x`$ be the unit vector with entry $`1`$ in position $`x`$, and $`0`$ elsewhere. Then it suffices to show that for all $`x,y`$

$`{\displaystyle \frac{1}{2d}}{\displaystyle \underset{a=1}{\overset{2d^2}{}}}{\displaystyle \underset{k=1}{\overset{d}{}}}x|E_{k,a}\rho E_{k,a}|y=\{\begin{array}{cc}1\hfill & \text{if }x=y\hfill \\ x|\rho |y\hfill & \text{if }xy\hfill \end{array}`$ (59)
The $`(a,k)^{\mathrm{th}}`$ term on the left side of (59) can be written as

$`x|E_{k,a}\rho E_{k,a}|y`$ $`=`$ $`x|\psi _{k,a}\psi _{k,a}|\rho |\psi _{k,a}\psi _{k,a}|y`$
$`=`$ $`{\displaystyle \underset{u,v=1}{\overset{d}{}}}x|\psi _{k,a}\psi _{k,a}|uu|\rho |vv|\psi _{k,a}\psi _{k,a}|y`$
Furthermore

$`x|\psi _{k,a}=x|G^kH^a|\theta ={\displaystyle \frac{1}{\sqrt{d}}}\mathrm{exp}\left[{\displaystyle \frac{2\pi ikx}{d}}\right]\mathrm{exp}\left[{\displaystyle \frac{2\pi iax^2}{2d^2}}\right]`$ (61)
Substituting into (3.1) gives

$`x|E_{k,a}\rho E_{k,a}|y`$ (62)
$`={\displaystyle \frac{1}{d^2}}{\displaystyle \underset{u,v=1}{\overset{d}{}}}u|\rho |v\mathrm{exp}\left[{\displaystyle \frac{2\pi ik}{d}}(x+vyu)\right]\mathrm{exp}\left[{\displaystyle \frac{2\pi ia}{2d^2}}(x^2+v^2y^2u^2)\right]`$
When the right side of (62) is substituted in (59), the sum over $`k`$ gives zero unless $`x+vyu`$ is an integer multiple of $`d`$. Since $`x,v,y,u`$ vary between $`1`$ and $`d`$, the only possible values are

$`x+vyu=0,d,d`$ (63)
Similarly, the sum over $`a`$ gives zero unless $`x^2+v^2y^2u^2`$ is an integer multiple of $`2d^2`$. In this case the only possibility is

$`x^2+v^2y^2u^2=0`$ (64)
Consider first the case that (63) gives

$`x+vyu=d`$ (65)
Let $`\gamma =y+u`$, then $`x+v=\gamma +d`$, and hence

$`2\gamma d`$ (66)
Elementary bounds then lead to

$`x^2+v^2>\gamma ^2>y^2+u^2`$ (67)
which shows that there can be no simultaneous solution of (64) and (65). A similar argument holds for the case

$`x+vyu=d`$ (68)
The remaining case in (63) can be written as

$`xu=yv`$ (69)
and also (64) can be written as

$`(xu)(x+u)=(yv)(y+v)`$ (70)
It follows that the only simultaneous solutions of the equations (69) and (70) are $`x=y,u=v`$ and $`x=u,y=v`$. Hence if $`xy`$ the left side of (59) gives $`x|\rho |y`$, while if $`x=y`$ the sum gives $`_uu|\rho |u=\mathrm{Tr}\rho =1`$. QED

Combining Lemma 10 and Lemma 11 we arrive at the convex decomposition of $`\mathrm{\Delta }_\lambda `$, which expresses the depolarizing channel in terms of the phase-damping channels $`\mathrm{\Phi }_\lambda ^{(a)}`$:

$`\mathrm{\Delta }_\lambda (\rho )`$ $`=`$ $`{\displaystyle \frac{\lambda }{1+(d1)\lambda }}{\displaystyle \frac{1}{2d}}{\displaystyle \underset{a=1}{\overset{2d^2}{}}}\mathrm{\Phi }_\lambda ^{(a)}(\rho )`$
$`+`$ $`{\displaystyle \frac{1\lambda }{1+(d1)\lambda }}{\displaystyle \frac{1}{2d^3}}{\displaystyle \underset{k=1}{\overset{d}{}}}{\displaystyle \underset{a=1}{\overset{2d^2}{}}}\left(G^{}\right)^k\mathrm{\Phi }_\lambda ^{(a)}(\rho )G^k`$
Furthermore each state $`|\psi _{k,a}`$ defined in (54) is uniform, and hence the phase-damping channels $`\mathrm{\Phi }_\lambda ^{(a)}`$ defined in (56) are also uniform. This completes the proof of Lemma 8.

### 3.2 $`d=2`$: the qubit depolarizing channel

It is useful to look in detail at the familiar case $`d=2`$. A general state $`\rho `$ can be written as a $`2\times 2`$ hermitian matrix

$`\rho =\left(\begin{array}{cc}a& c\\ \overline{c}& b\end{array}\right)`$ (72)
The depolarizing channel (1) acts by

$`\mathrm{\Delta }_\lambda (\rho )=\left(\begin{array}{cc}\lambda _+a+\lambda _{}b& \lambda c\\ \lambda \overline{c}& \lambda _{}a+\lambda _+b\end{array}\right),`$ (73)
where we have defined

$`\lambda _\pm ={\displaystyle \frac{1\pm \lambda }{2}}`$ (74)
Also the ‘intermediate’ channel $`\mathrm{\Omega }_\lambda `$ defined in (45) acts by

$`\mathrm{\Omega }_\lambda (\rho )=\left(\begin{array}{cc}\lambda _+a+\lambda _{}b& \lambda _+c\\ \lambda _+\overline{c}& \lambda _{}a+\lambda _+b\end{array}\right)`$ (75)
The first diagonal unitary matrix $`G`$ defined in (47) is just

$`G=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)=\sigma _z`$ (76)
So the first relation (49) becomes

$`\mathrm{\Delta }_\lambda ={\displaystyle \frac{2\lambda }{1+\lambda }}\mathrm{\Omega }_\lambda +{\displaystyle \frac{1\lambda }{1+\lambda }}{\displaystyle \frac{1}{2}}\left[\mathrm{\Omega }_\lambda +\sigma _z\mathrm{\Omega }_\lambda \sigma _z\right]`$ (77)
which can be easily verified using (73) and (75).

The second diagonal unitary matrix $`H`$ defined in (51) is now

$`H=\left(\begin{array}{cc}\mathrm{exp}\left[\pi i/4\right]& 0\\ 0& 1\end{array}\right)`$ (78)
There are eight phase-damping channels defined in (56). Four of these can be written in terms of the usual Pauli matrices:

$`\mathrm{\Phi }_\lambda ^{(2)}(\rho )=\mathrm{\Phi }_\lambda ^{(6)}(\rho )=\lambda _+\rho +\lambda _{}\sigma _y\rho \sigma _y`$ (79)
$`\mathrm{\Phi }_\lambda ^{(4)}(\rho )=\mathrm{\Phi }_\lambda ^{(8)}(\rho )=\lambda _+\rho +\lambda _{}\sigma _x\rho \sigma _x`$ (80)
The others can be written in terms of the following Pauli-type matrices:

$`\tau =\left(\begin{array}{cc}0& \mathrm{exp}\left[\pi i/4\right]\\ \mathrm{exp}\left[\pi i/4\right]& 0\end{array}\right),\overline{\tau }=\left(\begin{array}{cc}0& \mathrm{exp}\left[\pi i/4\right]\\ \mathrm{exp}\left[\pi i/4\right]& 0\end{array}\right)`$ (81)
The relations are

$`\mathrm{\Phi }_\lambda ^{(1)}(\rho )=\mathrm{\Phi }_\lambda ^{(5)}(\rho )=\lambda _+\rho +\lambda _{}\tau \rho \tau `$ (82)
$`\mathrm{\Phi }_\lambda ^{(3)}(\rho )=\mathrm{\Phi }_\lambda ^{(7)}(\rho )=\lambda _+\rho +\lambda _{}\overline{\tau }\rho \overline{\tau }`$ (83)
The second convex decomposition (57) now reads

$`\mathrm{\Omega }_\lambda ={\displaystyle \frac{1}{8}}{\displaystyle \underset{a=1}{\overset{8}{}}}\mathrm{\Phi }_\lambda ^{(a)}`$ (84)
There is a lot of redundancy in the final decomposition (3.1), which now has 24 terms on the right side. In fact $`\mathrm{\Omega }_\lambda `$ can be written as a convex combination of just two uniform phase-damping channels, namely

$`\mathrm{\Omega }_\lambda ={\displaystyle \frac{1}{2}}\left[\mathrm{\Phi }_\lambda ^{(2)}+\mathrm{\Phi }_\lambda ^{(4)}\right],`$ (85)
and this allows $`\mathrm{\Delta }_\lambda `$ to be written as a convex combination of just four phase-damping channels. There may be a similar redundancy in (3.1) for $`d>2`$.

## 4 The phase-damping channel

In this section we will establish Lemma 9 for the product channel $`\mathrm{\Phi }_\lambda I`$. Without loss of generality we will choose the basis $`=\{|i\}`$, so that $`\mathrm{\Phi }_\lambda `$ acts on a state by simply scaling all off-diagonal entries by the same factor $`\lambda `$, as in (29):

$`\mathrm{\Phi }_\lambda (\rho )=\lambda \rho +(1\lambda ){\displaystyle \underset{i=1}{\overset{d}{}}}|ii|i|\rho |i`$ (86)
The product channel $`\mathrm{\Phi }_\lambda I`$ acts on bipartite states $`\rho _{12}`$ defined on $`𝐂^d𝐂^d^{}`$ for some dimension $`d^{}`$. It will be convenient to view these states as $`d\times d`$ block matrices, where each block is itself a $`d^{}\times d^{}`$ matrix. Furthermore there is a convenient factorization of these blocks, which can be derived by the following observation. Let us write $`\sqrt{\rho _{12}}=(V_1\mathrm{}V_d)`$ where each $`V_i`$ is a $`dd^{}\times d^{}`$ matrix. Then we have

$`\rho _{12}`$ $`=`$ $`\left(\sqrt{\rho _{12}}\right)^{}\sqrt{\rho _{12}}`$
$`=`$ $`\left(\begin{array}{ccc}V_1^{}V_1& \mathrm{}& V_1^{}V_d\\ \mathrm{}& \mathrm{}& \mathrm{}\\ V_d^{}V_1& \mathrm{}& V_d^{}V_d\end{array}\right)`$
Recall the definition of $`\rho _2^{(i)}`$ in (33). With our choice of basis here, the matrix $`E_iI`$ is simply the orthogonal projector onto the $`i^{\mathrm{th}}`$ block on the main diagonal, hence

$`\rho _2^{(i)}=V_i^{}V_i`$ (88)
The key to deriving the bound (34) is to rewrite the factorization (4) as follows:

$`\rho _{12}`$ $`=`$ $`\left(\begin{array}{cccc}V_1^{}& 0& \mathrm{}& 0\\ 0& V_2^{}& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& V_d^{}\end{array}\right)M\left(\begin{array}{cccc}V_1& 0& \mathrm{}& 0\\ 0& V_2& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& V_d\end{array}\right)`$ (89)
where $`M`$ is the $`d\times d`$ block matrix

$`M=\left(\begin{array}{ccc}I^{}& \mathrm{}& I^{}\\ \mathrm{}& \mathrm{}& \mathrm{}\\ I^{}& \mathrm{}& I^{}\end{array}\right)`$ (90)
and $`I^{}`$ is the $`dd^{}\times dd^{}`$ identity matrix. Recall the state $`|\theta `$ defined in (53). Using this we can rewrite $`M`$ as the product state

$`M=d\left(|\theta \theta |\right)I^{}`$ (91)
Furthermore the simple action of the phase-damping channel $`\mathrm{\Phi }_\lambda `$ implies that it acts on (89) in the following way:

$`(\mathrm{\Phi }_\lambda I)(\rho _{12})`$ $`=`$ $`\left(\begin{array}{ccc}V_1^{}V_1& \mathrm{}& \lambda V_1^{}V_d\\ \mathrm{}& \mathrm{}& \mathrm{}\\ \lambda V_d^{}V_1& \mathrm{}& V_d^{}V_d\end{array}\right)`$
$`=`$ $`\left(\begin{array}{cccc}V_1^{}& 0& \mathrm{}& 0\\ 0& V_2^{}& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& V_d^{}\end{array}\right)(\mathrm{\Phi }_\lambda I)(M)\left(\begin{array}{cccc}V_1& 0& \mathrm{}& 0\\ 0& V_2& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& V_d\end{array}\right)`$
Let us define

$`A=\left(\begin{array}{cccc}V_1V_1^{}& 0& \mathrm{}& 0\\ 0& V_2V_2^{}& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& V_dV_d^{}\end{array}\right)`$ (94)
and

$`B=(\mathrm{\Phi }_\lambda I)(M)=d\left(\mathrm{\Phi }_\lambda (|\theta \theta |)\right)I^{}`$ (95)
Then $`(\mathrm{\Phi }_\lambda I)(\rho _{12})`$ has the same spectrum as the matrix $`A^{1/2}BA^{1/2}`$. Therefore

$`\mathrm{Tr}\left((\mathrm{\Phi }_\lambda I)(\rho _{12})\right)^p=\mathrm{Tr}\left(A^{1/2}BA^{1/2}\right)^p`$ (96)
Now we use the Lieb-Thirring inequality , which states that for all $`p1`$

$`\mathrm{Tr}\left(A^{1/2}BA^{1/2}\right)^p\mathrm{Tr}\left(A^{p/2}B^pA^{p/2}\right)=\mathrm{Tr}\left(A^pB^p\right)`$ (97)
The matrix $`A^p`$ is block diagonal:

$`A^p=\left(\begin{array}{cccc}(V_1V_1^{})^p& 0& \mathrm{}& 0\\ 0& (V_2V_2^{})^p& \mathrm{}& 0\\ \mathrm{}& & \mathrm{}& \mathrm{}\\ 0& & \mathrm{}& (V_dV_d^{})^p\end{array}\right)`$ (98)
Furthermore

$`B^p=d^p\left(\mathrm{\Phi }_\lambda (|\theta \theta |)\right)^pI^{}`$ (99)
Explicit calculation shows that the diagonal entries of $`d^p\left(\mathrm{\Phi }_\lambda (|\theta \theta |)\right)^p`$ are all equal to $`(1\lambda )^p+\left[(d\lambda +1\lambda )^p(1\lambda )^p\right]/d`$. Comparing this with (10), and substituting (98) and (99) into the right side of (97) we get

$`\mathrm{Tr}\left(A^pB^p\right)=d^{(p1)}\left(\nu _p(\mathrm{\Delta }_\lambda )\right)^p{\displaystyle \underset{i=1}{\overset{d}{}}}\mathrm{Tr}(V_iV_i^{})^p`$ (100)
Now recall (88), and also notice that for all $`i=1,\mathrm{},d`$

$`\mathrm{Tr}(V_iV_i^{})^p=\mathrm{Tr}(V_i^{}V_i)^p=\mathrm{Tr}\left(\rho _2^{(i)}\right)^p`$ (101)
Combining (97), (100) and (101) gives the bound (34). QED

## 5 The additivity of $`\chi ^{}`$

The proof of Theorem 2 uses the representation of $`\chi ^{}`$ as a min-max of relative entropy, combined with an entropy bound derived from Lemma 9. The relative entropy representation was derived by Ohya, Petz and Watanabe  and Schumacher and Westmoreland . Recall that the relative entropy of two states $`\rho `$ and $`\omega `$ is defined as

$`S(\rho ,\omega )=\mathrm{Tr}\rho (\mathrm{log}\rho \mathrm{log}\omega )`$ (102)
The OPWSW representation for the Holevo capacity of the channel $`\mathrm{\Psi }`$ is

$`\chi ^{}(\mathrm{\Psi })`$ $`=`$ $`\underset{\omega }{inf}\underset{\rho }{sup}S(\mathrm{\Psi }(\rho ),\mathrm{\Psi }(\omega ))`$ (103)
$`=`$ $`\underset{\rho }{sup}S(\mathrm{\Psi }(\rho ),\mathrm{\Psi }(\omega ^{}))`$ (104)
where the state $`\omega ^{}`$ that achieves the infimum in (103) is the optimal average input state from the channel. For the depolarizing channel this optimal average is $`(1/d)I`$, that is the totally mixed state. For the product channel $`\mathrm{\Delta }_\lambda \mathrm{\Psi }`$, the Holevo quantity $`\chi ^{}(\mathrm{\Delta }_\lambda \mathrm{\Psi })`$ can be upper bounded by choosing $`(1/d)I\omega ^{}`$ as the average input state. This leads to the following inequalities:

$`\chi ^{}(\mathrm{\Delta }_\lambda )+\chi ^{}(\mathrm{\Psi })\chi ^{}(\mathrm{\Delta }_\lambda \mathrm{\Psi })\underset{\tau _{12}}{sup}S((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))`$ (105)
In order to prove additivity we will combine this with the following result.

###### Lemma 12

For all bipartite states $`\tau _{12}`$,

$`S((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))\chi ^{}(\mathrm{\Delta }_\lambda )+\chi ^{}(\mathrm{\Psi })`$ (106)
Proof:

The left side of (106) can be rewritten as

$`S((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))`$ $`=`$ $`S\left((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})\right)`$
$`+`$ $`\mathrm{log}d\mathrm{Tr}\mathrm{\Psi }(\tau _2)\mathrm{log}\mathrm{\Psi }(\omega ^{})`$
where $`\tau _2`$ is the reduced density matrix of $`\tau _{12}`$. From here on we follow the steps in the proof of Theorem 3. First, by Lemma 7 we can assume without loss of generality that $`\tau _1=\mathrm{Tr}_2(\tau _{12})`$ is diagonal. Second, notice that the channel $`\mathrm{\Delta }_\lambda `$ appears on the right side of (5) only in the first term. Therefore Lemma 8 and concavity of the entropy imply that it is sufficient to establish the bound

$`S((\mathrm{\Phi }_\lambda \mathrm{\Psi })(\tau _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))\chi ^{}(\mathrm{\Delta }_\lambda )+\chi ^{}(\mathrm{\Psi })`$ (108)
where $`\mathrm{\Phi }_\lambda `$ is a uniform phase-damping channel and where $`\tau _1`$ is diagonal.

Next we apply (34) with $`\rho _{12}=(I\mathrm{\Psi })(\tau _{12})`$, and take the derivative at $`p=1`$ to get

$`S\left((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})\right)`$ $``$ $`S_{\mathrm{min}}(\mathrm{\Delta }_\lambda )\mathrm{log}d`$
$``$ $`{\displaystyle \underset{i=1}{\overset{d}{}}}x_i\mathrm{log}x_i+{\displaystyle \underset{i=1}{\overset{d}{}}}x_iS\left(\mathrm{\Psi }\left({\displaystyle \frac{1}{x_i}}\tau _2^{(i)}\right)\right)`$
where $`x_i=\mathrm{Tr}\left(\tau _2^{(i)}\right)`$, and as usual

$`\tau _2^{(i)}=\mathrm{Tr}_1\left[(E_iI)\tau _{12}\right]`$ (110)
Since $`\tau _1`$ is diagonal and $`\mathrm{\Phi }_\lambda `$ is uniform, it follows that

$`x_i={\displaystyle \frac{1}{d}}`$ (111)
for all $`i=1,\mathrm{},d`$, hence $`x_i\mathrm{log}x_i=\mathrm{log}d`$. Also recall the evaluation of $`\chi ^{}(\mathrm{\Delta }_\lambda )`$ in (11). Hence (5) can be written as

$`S\left((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\tau _{12})\right)\chi ^{}(\mathrm{\Delta }_\lambda )+\mathrm{log}d+{\displaystyle \frac{1}{d}}{\displaystyle \underset{i=1}{\overset{d}{}}}S\left(\mathrm{\Psi }\left(d\tau _2^{(i)}\right)\right)`$ (112)
Furthermore, since the projections $`E_i`$ in (110) constitute an orthonormal basis it follows that

$`{\displaystyle \underset{i=1}{\overset{d}{}}}\tau _2^{(i)}=\mathrm{Tr}_1\left[(II)\tau _{12}\right]=\tau _2`$ (113)
Therefore the left side of (108) can be rewritten as in (5) to get

$`S((\mathrm{\Phi }_\lambda \mathrm{\Psi })(\tau _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))`$ $`=`$ $`S\left((\mathrm{\Phi }_\lambda \mathrm{\Psi })(\tau _{12})\right)`$
$`+`$ $`\mathrm{log}d{\displaystyle \frac{1}{d}}{\displaystyle \underset{i=1}{\overset{d}{}}}\mathrm{Tr}\mathrm{\Psi }(d\tau _2^{(i)})\mathrm{log}\mathrm{\Psi }(\omega ^{})`$
Combining (5) with (112) we get

$`S((\mathrm{\Delta }_\lambda \mathrm{\Psi })(\rho _{12}),(1/d)I\mathrm{\Psi }(\omega ^{}))`$ $``$ $`\chi ^{}(\mathrm{\Delta }_\lambda )`$
$`+`$ $`{\displaystyle \frac{1}{d}}{\displaystyle \underset{i=1}{\overset{d}{}}}S(\mathrm{\Psi }(d\tau _2^{(i)}),\mathrm{\Psi }(\omega ^{}))`$
Recall that $`\mathrm{Tr}(d\tau _2^{(i)})=1`$. Therefore it follows from (103) that for each $`i=1,\mathrm{},d`$

$`S(\mathrm{\Psi }(d\tau _2^{(i)}),\mathrm{\Psi }(\omega ^{}))\chi ^{}(\mathrm{\Psi })`$ (116)
and hence (5) implies (108). QED

## 6 Conclusions and discussion

We have presented the proof of a long-conjectured property of the $`d`$-dimensional depolarizing channel $`\mathrm{\Delta }_\lambda `$, namely that its capacity for transmission of classical information can be achieved with product signal states and product measurements. This result follows as a consequence of several additivity results which we prove for the product channel $`\mathrm{\Delta }_\lambda \mathrm{\Psi }`$ where $`\mathrm{\Psi }`$ is an arbitrary channel. The principal result is the proof of the AHW conjecture for the matrix $`p`$-norm, for all $`p1`$, from which we deduce the additivity of minimal entropy and of the Holevo quantity. The argument presented here is a generalization of the method used earlier by the author to prove similar results for all unital qubit channels, and involves re-writing the depolarizing channel as a convex combination of other simpler channels, which we refer to as phase-damping channels.

If the additivity conjecture for the Holevo quantity is true for all channels, then there must be a general argument which can be used to provide a proof, and presumably this would give a different method of proof for Theorem 1. However it is known that the AHW conjecture is not true in general , and indeed it is an interesting problem to determine the class of channels for which it does hold. As a consequence, it may be that the method of this paper gives the most direct route to the proof of the AHW property for the depolarizing channel. It is expected that the same method can be applied to prove the AHW result for a class of $`d`$-dimensional channels, and this question is under study.

Acknowledgements The author is grateful to Wayne Sullivan and John Lewis for valuable discussions leading to the bound in Lemma 9. This work was partially supported by National Science Foundation Grant DMS–0101205, and by Science Foundation Ireland under the National Development Plan.


format: text

text:
# Numerical estimation of entropy loss on dimerization: improved prediction of the quaternary structure of the GCN4 leucine zipper

## I Introduction

Reduced models of a protein have been shown to provide a possible route for the estimation of the free energy of dimerization of relatively short coiled-coils . A key step in the calculation of the free energy of dimerization concerns the entropy loss upon bringing two monomer chains together to form the dimer. A practical method for the numerical estimation of this entropy loss by Monte Carlo simulation is discussed in this paper.

The focus of our work is on the calculation of free energies of dimerization, and in particular a re-analysis (using a subsequently improved model) of prior research about the folding thermodynamics of the GCN4 leucine zipper . Leucine zippers belong to the class of structural motifs that are known as coiled coils. Generically, they comprise right handed $`\alpha `$ helices wrapped around each other with a small left-handed super-helical twist . While leucine zippers can exist in both monomeric or dimeric form , GCN4 forms a dimer in the crystalline phase .

Numerical calculations of the free energy difference between the monomeric and dimeric forms of GCN4 have already been given in . In ref. , the free energy of the monomer was computed by transfer matrix methods, whereas the entropy change of dimerization was estimated by Monte Carlo methods. Two monomers were placed in a parallel configuration, and the configurational partition function was estimated by placing the two monomers in registry with one another, but considering different position of the relative starting point of each segment. The Entropy Sampling Monte Carlo method of Hao and Scheraga  was used in ref.  to sample the configurational space of both monomer and dimer. By using the values of the various terms that contribute to the entropy loss of dimerization given in , the equilibrium dimer fraction was computed as a function of temperature and monomer concentration. As a byproduct of the calculation, other equilibrium quantities such as the helical content of the monomer and dimer were obtained, as well as an analysis of the existence of possible folding intermediaries.

In our present work, we use the Replica Exchange Monte Carlo method  to obtain the canonical distribution of both monomer and dimer forms separately. Re-weighting methods are then used to calculate the respective free energies as a function of temperature . In the final step, a method is described to place both free energies on the same relative scale so that the free energy difference between the monomer and dimer forms can be computed. We believe that the method described in this paper allows a more accurate determination of the entropy loss of dimerization than previously reported, and therefore it allows a more accurate determination of the free energy difference between the monomeric and dimeric forms as a function of temperature and concentration.

## II The protein model

Our analysis is based on a a reduced model in which a protein is represented by a sequence of virtual bonds connecting effective particles . Each of these particles is assumed to be located at the center of mass of the side chain and backbone $`\alpha `$ carbon. The effective particles are embedded in a regular cubic lattice of fixed spacing that allows for a fairly accurate representation of the backbone of known protein structures. This geometric part of the model has been checked against all structures contained in the protein data bank . The observed root mean squared deviation between the lattice representation of any protein and its resolved structure is typically below 0.8 Å. The actual resolution of the model is of course lower, typically of the order of 2 Åfor small proteins.

Effective interactions (force fields) are introduced among the particles that include generic (sequence independent), and sequence specific contributions. The potentials associated with the generic type of interactions are defined so as to introduce a bias toward reasonable secondary structures. One such potential is introduced to account for the fact that proteins exhibit a characteristic bimodal distribution of neighbor residue distances, specially between the i-th and (i+4)-th residues. Configurations corresponding to the larger distance in the distribution are associated with proteins that exhibit either $`\beta `$-type or expanded coils, whereas the shorter distance corresponds to helices and turns . A second generic interaction further introduces a bias toward certain packing structures such as helices and $`\beta `$-type states. The first potential produces the required stiffness of the polypeptide chain, whereas the second provides for local cooperative motion during packing.

Sequence specific interactions are of three types, and include short range interactions, long range, pairwise interactions, and many body interactions. The short range pairwise potentials are of statistical origin and are fitted to nonhomologous reference structures. This is accomplished by considering the known distances between pairs of aminoacids that are separated by one through four bonds along the chain. Chirality is also introduced by considering an interaction between the i-th and (i+3)-th and (i+6)-th bonds to produce the correct pitch. Long range interactions are also of statistical origin, and are assumed to depend not only on relative distances between the effective atoms, but also on the relative orientation between the corresponding bonds (for example, residues of opposite charges are attractive when the corresponding bonds are parallel to each other, whereas the interaction is weak or repulsive when they are anti-parallel).

Finally, although multi-body interactions are implicitly included (as an unknown contribution) in the pairwise potentials determined from inter-residue distances and bond angles, two additional terms are added to model hydrophobic interactions and the known probability of a residue to have a given number of parallel and anti-parallel contacts. The hydrophobic potential is estimated from the surface exposure of a given side chain, i.e., of all possible contacts of a side chain, those that are not effectively occupied by contacts with neighboring chains. The second multi-body potential introduces a bias toward known propensities of various aminoacids to pack their side chains in parallel or anti-parallel orientation. Recent applications of the methodology include the improvement of threading based structure prediction , and direct ab-initio folding studies .

## III Monomer-dimer equilibrium and Monte Carlo method

At constant temperature and in thermodynamic equilibrium, the concentration of freely associating monomers and dissociating dimers is governed by the reaction constant, $`K=[D]/[M]^2`$, where $`[M]`$ and $`[D]`$ are the concentrations of monomers and dimers respectively. The mol fraction equilibrium constant can be expressed in terms of the canonical partition functions of both monomer and dimer as ,

$$K_x=c_0K=N\frac{Z_D}{Z_M^2}=N\frac{Q_D}{\sigma _DQ_M^2},$$
(1)
where $`N`$ is the total number of chains (i.e., two chains per dimer), $`c_0=N/V`$ the chain concentration in a system of volume $`V`$, and $`Z_M`$ and $`Z_D`$ the canonical partition functions of the monomer and dimer. The momenta degrees of freedom of both monomer and dimer can be integrated out from their partition functions, and the respective kinetic contributions cancel. We have therefore introduced the configurational partition functions $`Q_M`$ and $`Q_D`$, as well as the symmetry factor $`\sigma _D`$ that takes into account the indistinguishability of the two chains that form the dimer ($`\sigma _D=2!`$ in the present case) . The purpose of the present paper is an estimation of $`Q_M`$ and $`Q_D`$ by the Monte Carlo method.

Although it is in principle possible in to estimate the ratio of configurational partition functions in Eq. (1) by direct simulation involving coexisting monomers and dimers that transform into each other by association or dissociation, we have found that this is not practical. First, proteins are large molecules (even in the reduced representation used in our work and described in Appendix A), and only a small number of them can be placed in a computational cell. Second, the protein models employed lack very long range interactions, and hence there are large entropic contributions to the free energies which are difficult to sample accurately over the large configurational space of a protein. We instead follow the approach of refs.  which consists of two steps: an independent computation of the free energies of the monomer and dimer forms, followed by a transformation to a common reference state so that the free energy difference between the two can be estimated.

The configurational partition functions $`Q_M`$ and $`Q_D`$ are estimated by the Replica Exchange Monte Carlo method . Either a single monomer or a single dimer are placed in the computational cell, and a set of canonical Monte Carlo simulations are performed at set of prescribed neighboring temperatures. By conducting simulations involving only one monomer or one dimer we are explicitly assuming that at physiological concentrations monomer-dimer or dimer-dimer interactions can be neglected. The simulation to obtain $`Q_D`$ involves two identical monomers constrained during the course of the simulation to states in which there is at least one inter chain contact.

We consider a set of $`r`$ independent canonical simulations conducted in parallel at a set of neighboring inverse temperatures $`\beta _i,(i=1,\mathrm{},r)`$. At fixed intervals during the course of the simulation, two configurations at different temperatures (“replicas”) are chosen at random, and their respective temperatures exchanged with probability defined so as to preserve detailed balance as given by the canonical probability distribution. Additional details about the so-called Replica Exchange Monte Carlo method can be found in ref. . At each inverse temperature $`\beta _i`$, a sample of $`n_i`$ statistically independent configurations is collected, and the corresponding energy histograms $`h_i(E)`$ calculated with some arbitrary energy binning $`\mathrm{\Delta }E`$. Re-weighting is then used to estimate the partition function by simultaneously solving ,

$`p(E,\beta )`$ $`=`$ $`{\displaystyle \frac{\left(_{i=1}^rh_i(E)\right)e^{\beta E}}{_{i=1}^rn_ie^{\beta _iE+f_i}}},`$ (2)
$`e^{f_i}`$ $`=`$ $`\mathrm{\Delta }E{\displaystyle \underset{E}{}}p(E,\beta _i),`$ (3)
where $`f_i=\beta _iF(\beta _i)`$ is the dimensionless thermodynamic free energy at inverse temperature $`\beta _i`$. As is standard, free energy and canonical partition function are related through $`f_i=\mathrm{ln}Q(\beta _i)`$. Also, for latter reference, we note that the configurational density of states $`W(E)`$ is given by,

$$W(E)=p(E,\beta )e^{\beta E},$$

and therefore the configurational entropy can be obtained as,

$$\frac{S(E)}{k_B}=\mathrm{ln}\left(W(E)\mathrm{\Delta }E\right)$$

where $`k_B`$ is Boltzmann’s constant. Two sets of simulations are performed to yield monomer $`\{f_{M,i}\}`$ and dimer $`\{f_{D,i}\}`$ free energies at the set of inverse temperatures $`\{\beta _i\}`$. Note that each set is known up to an arbitrary additive constant.

In the remainder of this Section, we describe a number of transformations of the configurational partition functions $`Q_M`$ and $`Q_D`$ either for computational convenience or for the calculation of the proper reference state. Given the assumed integration of particle momenta, both $`Q_M`$ and $`Q_D`$ are expressed in terms of individual particle coordinates, and have dimensions of $`V^M`$ and $`V^{2M}`$ respectively, where $`M`$ is the number of aminoacids in the monomer. Elimination of rigid translation or rotation degrees of freedom from the configurational partition functions, for example, is usually accomplished by introducing the center of mass and principal axes of inertia of the molecule, and then relative coordinates for the individual particles. This requires either integrals over the corresponding canonical momenta, or the explicit consideration in the coordinate integrals of the appropriate transformation Jacobians. Since the Jacobians are configuration dependent, we follow instead earlier work , and conduct all our transformations on single particle coordinates alone thus obviating the need to introduce complicated Jacobian functions. For example, the elimination of the degrees of freedom associated with uniform translations is accomplished by eliminating the motion of particle 1 in one of the chains. The elimination of rigid rotation is partially accomplished by disallowing the rotation of the bond between particles 1 and 2 of one of the chains. Finally, in our estimate of entropy losses on dimerization (Section III A), we exclusively use phase space volumes of single particle coordinates to maintain the necessary dimensional consistency of Eq. (1) after all the transformations of both $`Q_M`$ and $`Q_D`$ that are described in this Section.

The configurational partition functions $`Q_M`$ and $`Q_D`$ are independent of the location of the molecule and of its orientation. The statistical accuracy of the simulation is greatly increased if those degrees of freedom that correspond to rigid translations and rotations are eliminated. The translational degree of freedom is eliminated by fixing the location of the first particle in the monomer chain, or of chain one in the dimer. Since the partition function is independent of this particle’s location, this coordinate can be integrated out to yield a factor of $`V`$ to the configurational partition function. The state of rigid rotation of the molecule can be specified by three angles, or the orientation of one axis plus a rotation around this axis. The orientation of the axis is defined by its azimuth $`\alpha (\pi ,\pi )`$ and a polar angle $`\theta (0,\pi )`$. The rotation around this axis is given by a third angle $`\gamma (\pi ,\pi )`$. Therefore the corresponding element of volume in configuration space is given by the triple integral $`𝑑\alpha d\mathrm{cos}(\beta )𝑑\gamma =8\pi ^2`$. This value can also be exactly factored out from the configurational partition function. However, given that the Monte Carlo method used employs multiple bond transitions, we have found it convenient to proceed somewhat differently. We disallow the Monte Carlo transition that corresponds to a two-bond change at the N-terminus, and therefore effectively eliminate the motion of particles 1 and 2 of the chain. While fixing the location of particle 1 still allows an exact calculation of the partition function, eliminating the motion of particle 2 introduces two approximations. The first one involves the factorization of the configuration space volume of particle 2. Since the motion of this particle is not independent of the motion of the rest of the chain, this factorization is only approximate. The degrees of freedom that are eliminated include a rotation around an axis defined by the bond vector between particles one and two, plus fluctuations in bond vector length. The elimination of the rotation is exact since the partition function is independent of the orientation of this axis, and yields a factor of $`𝑑\alpha d\mathrm{cos}(\beta )=4\pi `$ to the overall partition function. Factorization of the configuration space volume associated with bond length fluctuations is only approximate. We write,

$$QVV_2^{(1)}Q^{}V\frac{4\pi }{3}\left(R_{max}^3R_{min}^3\right)Q^{},$$
(4)
where $`Q^{}`$ is the partition function that is actually computed during the Monte Carlo simulation, and $`V_2^{(1)}`$ is the constant accessible volume for particle 2. The second approximation made involves the assumption that the second atom is free to move within a spherical shell centered in the first atom, of inner radius $`R_{min}=4.35`$ Å and outer radius $`R_{max}=7.94`$ Å. These two values are the smallest and largest bond distances allowed in the lattice model used. Finally, note that the computed partition function $`Q^{}`$ still contains a factor of $`2\pi `$ corresponding to the angle $`\gamma `$, the unrestricted rigid rotation of the molecule around the axis defined by the bond vector between particles one and two that is not eliminated during the simulation.

### A Reference state calculation

In order to introduce a common scale for the monomer and dimer free energy sets $`\{f_{M,i}\}`$ and $`\{f_{D,i}\}`$, we follow the method of ref. . At sufficiently high energies, one may assume that inter-chain interactions are negligible, and that the internal motions within each dimer chain are well approximated by those of the monomer. Therefore, and in this limit, the entropy of the dimer is approximately twice that of the monomer. This fact allows one to place the entropies from both monomer and dimer simulations in the same reference state at high energy, and hence to compute free energy differences between the two.

We briefly summarize here the steps taken for both monomer and dimer. The calculation of the internal entropy of the monomer is straightforward. In the dimer case, however, the contributions from the internal modes have to be separated from other degrees of freedom related to the relative position and orientation of the two chains. Since in the dimer simulation the two chains are constrained to have at least one contact, the computed entropies of the dimer at high energies still contain entropy losses due to this constraint that have to be estimated and subtracted to compute its internal entropy.

#### 1 Monomer

The monomer partition function $`Q_M^{}`$ computed by the Monte Carlo method still contains the contribution of one degree of freedom that is not associated with the internal motions of the particles, and that corresponds to a rigid rotation of the molecule around the axis defined by the bond between particles 1 and 2. From the Monte Carlo simulation and re-weighting, we obtain the probability density $`p^{}(E,\beta )`$ that corresponds to the partition function $`Q_M^{}`$ in Eq. (4). The corresponding entropy is,

$$\frac{S_M^{}(E)}{k_B}=\mathrm{ln}\left(p^{}(E,\beta )\mathrm{\Delta }E\right)+\beta E,$$

a quantity that is independent of the inverse temperature $`\beta `$. The internal entropy follows by subtracting the entropy of rotation of the azimuth angle of bond 2-3, and by adding the estimate given above for the radial part of the first bond,

$$\frac{S_{M,int}(E)}{k_B}=\mathrm{ln}\left(p^{}(E,\beta )\mathrm{\Delta }E\right)+\beta E\mathrm{ln}(2\pi )+\mathrm{ln}\frac{(R_{max}^3R_{min}^3)}{3}.$$
(5)
#### 2 Dimer

The dimer simulation is conducted by fixing the positions of atoms 1 and 2 of one of the chains. Therefore, a Monte Carlo estimate is obtained for $`Q_D^{}`$ as given in Eq. (4), As was the case for the monomer, we first define the entropy as estimated from the simulation by

$$\frac{S_D^{}(E)}{k_B}=\mathrm{ln}\left(p^{}(E,\beta )\mathrm{\Delta }E\right)+\beta E.$$

where $`p^{}(E,\beta )`$ corresponds to $`Q_D^{}`$ above. At sufficiently high energies, where the internal degrees of freedom of each chain are expected to become independent of the relative position and orientation of both chains, the total conformational density of states factors into a product involving the various contributions. In terms of the entropy, this factorization leads to the decomposition,

$`{\displaystyle \frac{S_{D,int}(E)}{k_B}}`$ $`=`$ $`\mathrm{ln}\left(p^{}(E,\beta )\mathrm{\Delta }E\right)+\beta E\mathrm{ln}(2\pi )+\mathrm{ln}{\displaystyle \frac{(R_{max}^3R_{min}^3)}{3}}\mathrm{ln}V_1^{(2)}(E)`$ (6)
$``$ $`\mathrm{ln}\left(\phi _\alpha (E)\phi _{\mathrm{cos}\beta }(E)\phi _\gamma (E)\right).`$ (7)
The quantity $`V_1^{(2)}(E)`$ is the accessible volume of particle 1 of chain 2 at energy $`E`$, and hence yields the accessible volume loss of dimerization. Its estimate during the Monte Carlo run is one of the main topics of this paper. Since the motion of the second chain relative to the first is constrained so that the number of inter-chain contacts is greater than zero, this volume will be in general much less than $`V`$. The quantity $`\phi _\alpha (E)`$ is the configuration volume available for the azimuth of bond 1 of chain 2, $`\phi _{\mathrm{cos}\beta }`$ for the cosine of the polar angle of bond 1 of chain 2, and $`\phi _\gamma `$ for the azimuth of bond 2 of chain 2. The product of the three represents the loss of rotational configuration space volume due to the formation of a dimer. If the second chain were to rotate freely relative to the first, we would have $`\phi _\alpha (E)\phi _{\mathrm{cos}\beta }(E)\phi _\gamma (E)=8\pi ^2`$. The value found is less that this upper bound, but it approaches $`8\pi ^2`$ as the energy of the dimer is increased. As was the case with $`V_2^{(2)}(E)`$, these three quantities also need to be estimated during the simulation.

#### 3 Configuration space volume estimation

The configuration space volumes $`V_1^{(2)}(E),\phi _\alpha (E),\phi _{\mathrm{cos}\beta }(E)`$ and $`\phi _\gamma (E)`$ have been estimated by using the method of coincidences . Let $`\mathrm{\Gamma }`$ be the volume of a certain region of configuration space. Consider a finite sample of configurations that are uniformly distributed in $`\mathrm{\Gamma }`$, and let $`\mathrm{\Gamma }_0\mathrm{\Gamma }`$ be a small coarse-graining volume in configuration space. The method involves computing the coincidence rate $`R`$ that a pair of configurations in the sample belongs to the same coarse-graining volume. If the configurations are uniformly distributed in $`\mathrm{\Gamma }`$, the probability of a coincidence is $`R=\mathrm{\Gamma }_0/\mathrm{\Gamma }`$. Therefore an estimate of $`R=n_c/n_t`$, where $`n_t`$ is the total number of pairs in the sample, and $`n_c`$ the total number of coincidences given $`\mathrm{\Gamma }_0`$, allows an estimation of $`\mathrm{\Gamma }`$.

In order to satisfy the conditions of the method, we first group all configurations (regardless of their temperature) according to their energy. Since all the configurations with the same energy are expected to occur with equal probability we calculate the coincidence rate $`R(E)`$ to estimate $`\mathrm{\Gamma }(E)`$ for the various magnitudes of interest ($`V_1^{(2)},\phi _\alpha ,\phi _{\mathrm{cos}\beta },`$ and $`\phi _\gamma `$).

We next note some limitations in the accuracy of the method. If $`\mathrm{\Gamma }_0`$ is not much smaller than $`\mathrm{\Gamma }`$, error is introduced as $`\mathrm{\Gamma }_0`$ will not generate a good covering set of $`\mathrm{\Gamma }`$, and it is likely that the method will overestimate the size of the region $`\mathrm{\Gamma }`$. On the other hand, if $`\mathrm{\Gamma }_0`$ is too small, the number of coincidences will be small, and the statistical error in the determination of $`R`$ is large. There is a third source of error associated with the sample size at each energy or, equivalently, the total number of pairs $`n_t(E)`$ . The number of coincidences can be estimated as,

$$n_c(E)\frac{1}{2}\left(\frac{n_t(E)}{k}\right)^2\frac{\mathrm{\Gamma }_0}{\mathrm{\Gamma }},$$

with $`\mathrm{\Gamma }_0`$ such that all $`n_t(E)`$ configurations have been distributed among $`k`$ groups. Therefore,

$$\mathrm{\Gamma }\frac{1}{2}\left(\frac{n_t(E)}{k}\right)^2\frac{\mathrm{\Gamma }_0}{n_c(E)},$$

so that for a fixed minimum $`n_c(E)`$ to insure adequate statistics of the coincidence rate, the estimated value of $`\mathrm{\Gamma }`$ is bounded by $`n_t^2(E)`$. Therefore sufficiently large samples are needed at each energy if the corresponding value of $`\mathrm{\Gamma }`$ is large. We will further illustrate these limitations in Section IV.

#### 4 Reference entropy difference

In order to place both the monomer and dimer in the same reference state, we require that in the limit of high $`E`$,

$$S_{D,int}(E)+S_0=2S_{M,int}(E/2),$$
(8)
where $`S_0`$ is a constant, independent of $`E`$, $`S_{D,int}(E)`$ is given by Eq. (6) and $`S_{M,int}(E)`$ by Eq. (5). The quantity $`S_0`$ is determined numerically as shown in Section IV.

Once the constant $`S_0`$ has been determined, the free energy and partition function of the dimer are re-scaled according to

$$f_i^{\prime \prime }=f_i^{}S_0,Q_D^{\prime \prime }=Q_D^{}e^{S_0}.$$
(9)
We can now compute the equilibrium constant $`K_x`$ by substituting Eq. (4) for both monomer and dimer into Eq. (1), but using the rescaled dimer partition function $`Q_D^{\prime \prime }`$ defined in Eq. (9) instead of $`Q_D^{}`$,

$$K_x=\frac{N}{VV_2^{(1)}}\frac{Q_D^{\prime \prime }}{\sigma _DQ_M^{\mathrm{\hspace{0.33em}2}}}=\frac{N}{VV_2^{(1)}}\frac{e^{S_0}Q_D^{}}{\sigma _DQ_M^{\mathrm{\hspace{0.33em}2}}}.$$
(10)
With this re-definition of the dimer partition function, both $`Q_D^{\prime \prime }`$ and $`Q_M^{}`$ are referred to the same reference state, and hence absolute values of $`K_x`$ can be given.

In terms of the free energies of the monomer and dimer that are obtained from the Monte Carlo calculation after re-weighting (Eqs. (2) and (3)), Eq. (10) leads to,

$$\mathrm{ln}K_x(\beta _i)=\mathrm{ln}c_0+S_0f_{D,i}^{}\mathrm{ln}\sigma _D\mathrm{ln}V_2^{(1)}+2f_{M,i}^{}.$$
(11)
If $`\mathrm{ln}K_x>0`$ the dimer is prevalent.

## IV Results

The method described in Section III has been tested on the GCN4 leucine zipper (a 31 residue segment with the characteristic heptad repeat sequence of leucine zippers). The oligomerization equilibrium of the wild type has been addressed both experimentally  and computationally , as well as that of several of its mutant forms . Due to the short length of the sequence, and the simplicity of its secondary structure, numerous computational studies have addressed various aspects of the oligomerization process in GCN4, including dimer and multi-mer equilibria , the stability of several of its sub-domains , oligomeric equilibrium of several of its mutant forms , and other parameters of the coiled coil such as the helical content as a function of temperature and a van’t Hoff enthalpy analysis to reveal the adequacy of a two state assumption for the dimerization process .

We have extended the analysis of  in two directions. First we use a Replica Exchange Monte Carlo method instead of the Entropy Sampling Monte Carlo method of that reference as the former provides a faster rate of convergence to the equilibrium distribution of the dimer form. Second, we extend the method of calculation of the various entropy losses upon dimerization, and show their strong dependence on the energy of the configuration, a dependence that was not taken into account in previous studies.

The results shown are based on two long runs for the monomer and dimer forms respectively. Initial configurations were chosen close to the native state, but first equilibrated at constant temperature. Several runs with different initial conditions yielded essentially identical results for the various thermodynamic quantities presented, although none of the dimer simulations involved an initial condition in a manifestly anti-parallel configuration. The monomer runs involved $`2\times 10^6`$ independent configurations or steps after equilibration, with one replica exchange attempted every 500 steps. Quantities for analysis were collected every 250 steps. The dimer simulation comprises two identical, and initially parallel chains with at least one contact between them . The simulation in this case is conducted by rejecting all bond moves that would result in no contacts between the chains. The run for the dimer involved $`1.3\times 10^6`$ configurations, with the same frequency of analysis and of replica exchange. In both cases, twenty independent replicas were run in parallel at dimensionless temperatures in the range $`T=0.51.45`$ in increments of 0.05. In the low temperature range, the root mean squared deviation (RMSD) between the estimated location of the $`\alpha `$ carbons in the model and the native configuration is of the order of 3 Å(see Fig. 1). We note that this RMSD range constitutes a prediction, and was not enforced during the course of the simulation.

The analysis presented is based on energy histograms, and the subsequent re-weighting described in Section III. The entire range of energies sampled by the monomer and the dimer during the course of the simulation was divided into 100 equal bins, so that in dimensionless units $`\mathrm{\Delta }E0.2926`$ for the monomer and $`\mathrm{\Delta }E0.3767`$ for the dimer.

In the case of the dimer, the location of all individual particles was also recorded every 250 steps in order to estimate the configuration space volumes $`V_1^{(2)},\phi _\alpha ,\phi _{\mathrm{cos}\beta }`$ and $`\phi _\gamma `$. The results presented for $`V_1^{(2)}`$ are based on the spatial coordinates of particle 15 of chain 2 (the chain that is free to move within the computational cell). Substantially identical results follows from an analysis of any other particles in the chain, except for those in the immediate vicinity of the N- or C- termini.

The configuration space volume $`\phi _\alpha `$ is estimated from the azimuth of the bond between particles 15 and 16 of chain 2, and $`\phi _{\mathrm{cos}\beta }`$ follow from $`\mathrm{cos}\beta `$, $`\beta `$ being the polar angle of this bond. Finally, $`\phi _\gamma `$ is obtained from the azimuth distribution of the bond between particles 15 and 16. In a freely rotating molecule, $`\alpha `$ is uniformly distributed in $`(\pi ,\pi )`$, $`\mathrm{cos}\beta `$ in $`(1,1)`$, and $`\gamma `$ in $`(\pi ,\pi )`$, resulting in a combined conformational space volume for rigid rotation of $`\phi _\alpha \phi _{\mathrm{cos}\beta }\phi _\gamma =8\pi ^2`$.

Figure 2 shows our results for $`V_1^{(2)}`$ with the same energy bin size $`\mathrm{\Delta }E`$ used to construct the histogram. The coarse-graining volume $`\mathrm{\Gamma }_0=\mathrm{\Delta }x\mathrm{\Delta }y\mathrm{\Delta }z`$ has been obtained by defining $`\mathrm{\Delta }x=(x_{max}x_{min})/\delta `$, and similarly for $`\mathrm{\Delta }y`$ and $`\mathrm{\Delta }z`$. $`x_{min}`$ and $`x_{max}`$ are the smallest and largest values of $`x_{15}`$, the $`x`$ coordinate of particle 15, for each particular energy bin. We present our results for a range of values of $`\delta `$ in Fig. 2. If $`\delta `$ is too large, the coarse-graining volume is small, and the number of configurations for a given energy $`n_t(E)`$ is also small. As discussed in Section III, this leads to underestimate the accessible volume. The value of $`V_1^{(2)}`$ is seen to increase with decreasing $`\delta `$, becomes approximately independent of $`\delta `$ in some range, and then further increases with decreasing $`\delta `$. If $`\delta `$ is too small, the shape of the region being sampled cannot be accurately reproduced with this coarse $`\mathrm{\Gamma }_0`$. Note that the values of $`V_1^{(2)}`$ at low energies are the most difficult to estimate, presumably because the shape of the region in configuration space is not as smooth as that at higher energies. However since the procedure leading to the computation of the reference entropy relies only on the region of high energies, this inaccuracy does not represent a significant limitation to our results.

The behavior just described is qualitatively similar to that shown in Fig. 3 corresponding to the rotation volume $`\phi _\alpha \phi _{\mathrm{cos}\beta }\phi _\gamma `$. In this case we define $`\mathrm{\Gamma }_0=\mathrm{\Delta }\alpha \mathrm{\Delta }(\mathrm{cos}\beta )\mathrm{\Delta }\gamma `$ with the same definition of $`\mathrm{\Delta }\alpha ,\mathrm{\Delta }(\mathrm{cos}\beta )`$ and $`\mathrm{\Delta }\gamma `$ in terms of the quantity $`\delta `$. The figure also shows (solid line) the value $`8\pi ^2`$ that corresponds to free rotation of chain 2 relative to chain 1. As can be seen from the figure, the values obtained approach this limit at high energies.

The constant $`S_0`$ of Eq. (8) required to place both the monomer and dimer free energies in the same scale is obtained directly from Eqs. (5) and (6), as shown in Fig. 4. In this figure we plot $`S_{D,int}(E)+S_0`$ and $`2S_{M,int}(E/2)`$ with $`S_0`$ adjusted graphically so that the two curves coincide at large $`E`$. Note that both curves superimpose to a good accuracy for a range of energies, indicating the consistency of the approach.

We show next our results for $`K_x`$ in Fig. 5, with $`K_x`$ defined in Eq. (10), as a function of the dimensionless temperature. When $`K_x=1`$ the mol fractions of the monomer and dimer forms are equal ($`x_M=x_D=0.5`$). At low temperatures $`K_x>1`$, indicating a prevalence of the dimer form, and the reverse is true at high temperatures. For the sake of illustration, the figure shows the values of $`K_x`$ at two different concentrations $`c_0=10\mu M`$ and $`c_0=1mM`$.

It is also intersting to examine the contact map of the dimer phase as given by the simulation. As discussed above, the only constraint in the simulation is that there be at least one contact between the two chains. Therefore the question arises as to whether the dimer retains a significant fraction of native contacts in the vicinity of the transition temperature, or whether there is a significant fraction of out of register dimer configurations that are structurally very different from the native state. In order to answer this question, we proceed as follows: A contact between residues belonging to different chains is considered native if it appears in the contact map of the native protein in its lattice representation. With this definition, GCN4 has 10 inter chain native contacts. We then calculate the ensemble average of the fraction of configurations that have at least 50 % native contacts. The results are shown in Fig. 6 as a function of temperature. This fraction approaches one at low temperature, changes quickly around the transition region, decaying to zero at high tempertures. We also shown in this figure the equilibrium mol fraction of the dimer form for the same two concentrations shown in Fig. 5. From the figure, we conclude that in this range of physiological concentrations, the decrease in the fraction of native contacts as given by our model can be mainly attributed to the appearance of the monomer form, and not to a significant contribution from out of register dimers.

To conclude, our results show that it is possible to calculate the entropy loss of dimerization corresponding to the GCN4 monomer-dimer equilibrium without any restrictions to the motion of the individual chains. Previous research on this system differed from ours in that the entropy loss was estimated by restricting the conformation space of the dimer, thus resulting in low values of the entropy (compare the value of $`V_1^{(2)}=67.6`$Å<sup>3</sup> given in Table 1 of ref. , and the values shown in Fig. 2). Despite the fact that we have allowed sampling of the full conformational space of the dimer, the results obtained confirm that it is possible to obtain the free energy of dimerization of GCN4 at physiological concentrations by using a reduced model of the protein and Monte Carlo simulations.

## Acknowledgments

This research has been supported by the National Science Foundation, grant No. 9986019. JV is also supported by the National Institutes of General Medical Sciences, grant No. GM64150-01.

## A Lattice protein model and interaction force parameters

The model protein used in this work employs a reduced representation of the protein backbone on a regular lattice. The model comprises a sequence of bonds connecting particles located at the center of mass of the corresponding residue and backbone $`\alpha `$ carbon. The particles are then placed in a three dimensional simple cubic lattice with spacing of 1.45 Å. Further details on this model can be found in ref.

A sequence of configurations is generated by a Monte Carlo scheme with Metropolis updating. The method employs three different types of individual transitions. In the first case, a single particle and its two corresponding bonds are selected for an attempted update. The second type of transition involves three consecutive bonds and the corresponding two adjacent particles. The third involves a rigid translation of a small fragment of the chain comprising three particles, and the ensuing rearrangement of the end bonds. These transitions are attempted sequentially for all the bonds in the chain. Two separate transitions are also included to adjust the position of the N- and C-termini particles. The set of all these attempted transitions constitutes a Monte Carlo step (MCS). Further details about the transitions used in the Monte Carlo updating can be found in ref. .

Interaction forces can be grouped into generic and sequence specific. The former are sequence independent and lead to protein-like packing, whereas the latter are derived from a statistical analysis of the protein database, and explicitly depend on the identity of the aminoacids involved. We next list the values of the various parameters used in our calculations. Sequence dependent short range interactions are defined by Eq. (1) of . We use a common multiplicative factor in our calculations $`\epsilon _{short}=0.325`$ (this factor is explicitly shown in Eq. (12) of  with a value of 0.75 instead). A three-body potential that is sequence specific is also used, with an amplitude $`\epsilon _{3b}=0.25`$. The generic, short range conformational biases of Eqs. (3), (4), and (5) of  involve $`\epsilon _{gen}=1.25`$. Hydrogen bonding energies within the main chain are also included, with an amplitude $`\epsilon _{Hbond}=0.325`$ in Eq. (7) of . Long range and sequence dependent interactions are modeled by a set of square well potentials as described in Eq. (9) of . We have chosen $`E^{rep}=4`$ in Eq. (8) of that reference, and a common multiplicative factor $`\epsilon _{pair}=2.0`$ (to be compared with the factor of 1.25 in Eq. (12) of ). Two additional multibody potentials are introduced to include hydrophobic effects, and preferences for parallel or anti-parallel packing among the residues. We define as the scale of Eq. (10) in  $`\epsilon _{surface}=0.75`$ (instead of the value 0.5 shown in Eq. (12) of that reference). Finally, we have used a factor $`\epsilon _{multi}=0.75`$ in Eq. (11) of  (instead of the value 0.5 in Eq. (12)).


text:
# References

Bosonic stimulation and the irreproducibility of condensate fragmentation

N. Kumar

Raman Research Institute, Bangalore 560080, India

## Abstract

It is pointed out that the quantum statistical phenomenon of Bosonic stimulation, inherent to Bose statistics and the associated Bose-Einstein correlation, can be effectively mapped on to the statistical problem of the Polya urn scheme. Thus, we predict an irreproducibility for the limiting non-degenerate values of the relative populations of two, or more equivalent but separated condensates resulting from the fragmentation of a given source condensate. Experiments are proposed that should verify this prediction $``$ using the dilute gas Bose-Einstein condensates, or a pulsed degenerate multi-mode laser, where one would look for the run-to-run (or the pulse-to-pulse) fluctuations of the relative populations.

It is a deep, but rather well known result that the indistinguishability of identical Bosons, and the associated quantum (Bose) statistics, implies a certain non-classical Bose-Einstein correlation.1,2 Thus, the indistinguishable Bose particles emitted from a given source, or from proximate sources, are so correlated in space-time, or momentum-energy, as to imply enhanced probability of emission of the particles with small relative momentum. A direct manifestation of this Bosonic stimulation,1,3 wherein a Bose particle is preferentially scattered (stimulated scattering) into the final single particle state which has a higher pre-existing Bosonic population. This derives straightforwardly from the usual Bose factor ($`n`$+1) in $`a^+|n=\sqrt{n+1}|n+1`$. It is now easy to see that, from the point of view of the relative population statistics, this phenomenon of Bosonic stimulation maps on to that associated with the Polya urn scheme.4 Recalling briefly, here one begins with one black and one white ball in the urn. The rule of the game then is to draw a ball at random from the urn, note its colour, and then replace the ball so drawn back into the urn, but with one additional ball of the same colour. The relative populations of the black and the white balls fluctuate initially, but eventually iterate away to a limiting value for the given run deterministically $``$ i.e., with the fluctuations about the limiting value staying within arbitrarily small tolerance chosen as $`n\mathrm{}`$, where $`n`$ numbers the draws. The limiting value, however, varies randomly from one run to another $``$ hence the irreproducibility of the limiting non-degenerate value. It is the vagaries of the initial few draws that determine the limiting value for a given run. Now, it is clear from the above that the rule to replace the ball so drawn back into the urn, but with one additional ball of the same colour simulates the ($`n`$+1) factor of the Bosonic stimulation. Its relation to a fragmentation of a BEC is now clear.

More precisely, consider a Bose-Einstein condensate as a source from which the Bosons could be coupled out into two equivalent but well separated traps. Bosonic stimulation then determines the occupation numbers of the two traps, as given by

$$P_{t+1}(n_1,n_2)=P_t(n_11,n_2)\left(\frac{n_1}{n_1+n_2+1}\right)+P_t(n_1,n_21)\left(\frac{n_2}{n_1+n_2+1}\right),$$
(1)
where $`P_t(n_1,n_2)`$ is the probability of occupation numbers ($`n_1,n_2`$) for the two traps (1,2) at time t. The time here is discretized so that one particle is out-coupled in the unit time step.

This is the probabilistic description of the Polya urn scheme. Clearly, as $`t\mathrm{}`$, one settles down to relative populations distributed uniformly in the interval zero to one. The result readily generalizes to more than two traps $``$ multifragmentation.

Fragmented B-E condensate is normally understood in terms of repulsive inter-particle interaction. As shown by Nozieres3, for a repulsive inter-particle interaction the Bosonic exchange makes a fragmentation of the condensate energetically costly in equilibrium. What we are discussing here, however, is the purely quantum statistical effect of Bosonic stimulation on the process of fragmentation $``$ a partitioning into the two traps which do not interact dynamically. The non-degenerate irreproducible relative populations are a type of partitioning noise $``$ but from one run to another.

Yet another experimental possibility is to consider a degenerate, two-mode laser which is pulsed. Here again, the photon released from the source (a stimulated emission from inverted atomic population) can populate the two degenerate cavity modes, from which the photons are coupled out eventually. Exactly the same considerations apply, with the only difference that in a pulsed laser the pulse-to-pulse fluctuations give directly the statistics of the limiting non-degenerate values of the relative populations. That should make it experimentally easier to realize.

Finally, we note that our work addresses an aspect of the Bose-Einstein correlation that enters many other fields $``$ including astronomy, the celebrated Hanbury-Brown and Twiss2,5 (HBT) effect. The dilute gaseous BEC offers new possibilities here.

The author would like to thank Hema Ramachandran and Andal Narayanan for discussion on possible experimental demonstration of the above effect. It is a pleasure to thank Robert M. May for a discussion on the Polya urn scheme.


text:
# Critical Temperature of an Interacting Bose Gas in a Generic Power-Law Potential

## References

1. M.H. Anderson, J.R. Ensher, M.R. Matthews, C.E. Wieman, and E.A. Cornell, Science 269 (1995) 189.

2. K.B. Davis, M.O. Mewes, M.R. Andrews, N.J. van Druten, D.S. Drufee, D.M. Kurn, and W. Ketterle, Phys. Rev. Lett. 75 (1995) 3969.

3. C.C. Bradley, C.A. Sackett, J.J. Tollet, and R.G. Hulet, Phys. Rev. Lett. 75 (1995) 1687.

4. F. Dalfovo, S. Giorgini, L.P. Pitaevskii, and S. Stringari, Rev. Mod. Phys. 71 (1999) 463.

5. V. Bagnato, D.E. Pritchard and D. Kleppner, Phys. Rev. A 35 (1987) 4354.

6. S. Grossmann and M. Holthaus, Phys. Lett. A 208 (1995) 108.

7. S. Giorgini, L. Pitaevskii and S. Stringari, Phys. Rev. A 54 (1996) 4633.

8. B. Pozzi, L. Salasnich, A. Parola and L. Reatto, Eur. Phys. Jour. D 11 (2000) 367.

9. W. Ketterle and D.E. Pritchard, Phys. Rev. A 46(1992) 4051.

10. P.W.H. Pinkse, A. Mosk, M. Weidenmuller, M.W. Reynolds, T.W. Hijmans, and J.T.M. Walraven, Phys. Rev. Lett. 78 (1997) 990.

11. D.M. Stamper-Kurn, H.J. Miesner, A.P. Chikkatur, S. Inouye, J. Stenger, and W. Ketterle, Phys. Rev. Lett. 81 (1998) 2194.

12. L. Salasnich, J. Math. Phys. 41 (2000) 8016.

|  | $`1`$ | $`2`$ | $`3`$ | $`4`$ | $`5`$ | $`6`$ | $`7`$ | $`8`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| $`c(n)`$ | $`0.56`$ | $`0.94`$ | $`1.14`$ | $`1.25`$ | $`1.31`$ | $`1.35`$ | $`1.04`$ | $`1.00`$ |
| $`d(n)`$ | $`0.64`$ | $`1.33`$ | $`0.23`$ | $`2.41`$ | $`7.24`$ | $`\mathrm{}`$ | $`41.57`$ | $`42.96`$ |
| $`e(n)`$ | $`0.40`$ | $`0.50`$ | $`0.29`$ | $`1.11`$ | $`2.45`$ | $`\mathrm{}`$ | $`2.96`$ | $`1.61`$ |
| $`\frac{2n}{3(n+2)}`$ | $`2/9`$ | $`1/3`$ | $`2/5`$ | $`4/9`$ | $`10/21`$ | $`1/2`$ | $`14/27`$ | $`8/15`$ |
| $`\frac{n}{3(n+2)}`$ | $`1/9`$ | $`1/6`$ | $`1/5`$ | $`2/9`$ | $`5/21`$ | $`1/4`$ | $`7/27`$ | $`4/15`$ |

TAB. 1. Numerical coefficients as a function of the power-law exponent $`n`$. In the limit $`n\mathrm{}`$ one gets: $`c(n)1.27`$, $`d(n)0`$, $`e(n)0.29`$, $`2n/(3(n+2))2/3`$, and $`n/(3(n+2))1/3`$.