format: text

text:
# A Network-Flow Technique for Finding Low-Weight Bounded-Degree Spanning Trees11footnote 1A preliminary version of this paper appeared in the Proceedings of the 5th International Integer Programming and Combinatorial Optimization Conference (IPCO), June 1996, pages 105-117.

## 1 Introduction

Given a complete graph with edge weights satisfying the triangle inequality, and a degree bound for each vertex, we consider the problem of computing a low-weight spanning tree in which the degree of each vertex is at most its given bound. In general, it is NP-hard to find such a tree. There are various practical motivations: the problem arises in the context of VLSI layout and network design  (such as in the Bellcore software FIBER OPTIONS, used for designing survivable optimal fiber networks). The special case of only one vertex with a degree-constraint has been examined ; a polynomial time algorithm for the case of a fixed number of nodes with a constrained degree was given by Brezovec et al. . Computational results for some heuristics for the general problem are presented in . Papadimitriou and Vazirani  raised the problem of finding the complexity of computing a minimum-weight degree-4 spanning tree of points in the plane. Some geometric aspects are considered in .

In this paper, we consider modifying a given spanning tree $`T`$, to meet the degree constraints without increasing its weight considerably. We introduce a novel network-flow based algorithm that does this optimally in the following sense: if for some algorithm a worst-case performance guarantee can be proved that is solely a function of the topology and edge weights of $`T`$, then that performance guarantee also holds for our algorithm. We prove this by showing that our algorithm finds the optimal solution for graphs in which the weight of each edge $`(u,v)`$ equals the cost of the $`uv`$ path in $`T`$.

We also show the following more concrete performance guarantee: If the degree constraint $`d(v)`$ for each $`v`$ is at least $`2`$, our algorithm finds a tree whose weight is at most the weight of $`T`$ times

$$2\mathrm{min}\left\{\frac{d(v)2}{deg_T(v)2}:deg_T(v)>2\right\},$$

where $`deg_T(v)`$ is the initial degree of $`v`$. For instance, the degree of each vertex $`v`$ can be reduced by nearly half, to $`1+deg_T(v)/2`$, without increasing the weight of the tree by more than $`50\%`$. (For comparison, note that a factor of 2 is straightforward with standard shortcutting techniques.) We also describe linear-time algorithms that achieve this ratio.

This performance guarantee is optimal in the sense that for any $`Dd2`$, if $`T`$ is a complete rooted $`(D1)`$-ary tree with unit edge weights and the edge weights in $`G`$ are those induced by paths in $`T`$, then the weight of any spanning tree with maximum degree $`d`$ is at least the weight of $`T`$ times $`2\frac{d2}{D2}o(1)`$.

The restriction $`d(v)2`$ is necessary to obtain constant performance bounds. Consider the case when $`T`$ is a simple path of unit weight edges, with the remaining edge weights again induced by $`T`$. Any spanning tree in which all but one vertex has degree one is heavier than $`T`$ by a factor of $`\mathrm{\Omega }(n)`$, the number of vertices in $`T`$.

For many metric spaces, graphs induced by points in the space have minimum spanning trees of bounded maximum degree. In such cases our algorithms can be used to find spanning trees of even smaller degree with weight bounded by a factor strictly smaller than 2 times the weight of a minimum spanning tree (MST). For example, in the $`L_1`$ metric, a degree-4 MST can be found , so that we can find a degree-3 tree with weight at most $`1.5`$ times the weight of an MST. We discuss similar results for the $`L_1`$, $`L_2`$, and $`L_{\mathrm{}}`$ norms. For some of these norms, this improves the best current performance guarantees.

Finally, we disprove the following conjecture of : “In Euclidean graphs, perhaps a Traveling Salesman path of weight at most $`(2\epsilon )`$ times the minimum spanning-tree weight always exists…”

Our algorithms modify the given tree by performing a sequence of adoptions. Our polynomial-time algorithm performs an optimal sequence of adoptions. Adoptions have been previously used to obtain bounded-degree trees in weighted graphs . The main contributions of this paper are a careful analysis of the power of adoptions and a network-flow technique for selecting an optimal sequence of adoptions. The method yields a stronger performance guarantee and may yield better results in practice. The analysis of adoptions shows that different techniques will be necessary if better bounds are to be obtained.

In the full version of their paper, Ravi et al. \[16, Thm. 1.9\] (if slightly generalized and improved222To obtain the improved bound one has to change the proof slightly by upper bounding $`c(v_1v_2)c(vv_2)`$ by $`c(vv_1)`$ and not $`c(vv_2)`$ as is done in .) gave an algorithm with a performance guarantee of

$$2\mathrm{min}\left\{\frac{d(v)2}{deg_T(v)1}:vV,deg_T(v)>2\right\}$$

provided each $`d(v)3`$. The performance guarantee of our algorithm is better.

In Euclidean graphs (induced by points in $`\mathrm{IR}^d`$), minimum spanning trees are known to have bounded degree. For such graphs, Khuller, Raghavachari and Young  gave a linear-time algorithm to find a degree-3 spanning tree of weight at most $`5/3`$ times the weight of a minimum spanning tree. For points in the plane, the performance guarantee of their algorithm improves to $`1.5`$; if the tree is allowed to have degree four, the ratio improves further to $`1.25`$.

In unweighted graphs, Fürer and Raghavachari  gave a polynomial-time algorithm to find a spanning tree of maximum degree exceeding the minimum possible by at most one. In arbitrary weighted graphs, Fischer  showed that a minimum spanning tree with maximum degree $`O(\delta ^{}+\mathrm{log}n)`$ can be computed in polynomial time, where $`\delta ^{}`$ is the minimum maximum degree of any minimum spanning tree. He also provided an algorithm that finds a minimum spanning tree with degree $`k(\delta ^{}+1)`$ where $`k`$ is the number of distinct edge weights.

## 2 Adoption

Fix the graph $`G=(V,V\times V)`$ and the edge weights $`w:V\times V\mathrm{IR}`$. The algorithm starts with a given tree $`T`$ and modifies it by performing a sequence of adoptions. The adoption operation (illustrated in Figure 1) is as follows:

| $`\text{Adopt}(u,v)`$ |
| --- |
|  | Precondition: Vertex $`v`$ has degree at least two in the current tree. |
| 1 | Choose a neighbor $`x`$ of $`v`$ in the current tree other than the neighbor on the |
|  | current $`uv`$ path. |
| 2 | Modify the current tree by replacing edge $`(v,x)`$ by $`(u,x)`$. |

The effect of $`\text{Adopt}(u,v)`$ is that $`u`$ adopts a neighbor of $`v`$. For the purpose of the algorithm, it makes no difference which neighbor of $`v`$ is adopted, since the differences in cost are hidden by assigning $`w(u,v)`$ as the cost of the adoption as shown below. $`\text{Adopt}(u,v)`$ decreases the degree of $`v`$ by one, at the expense of increasing the degree of $`u`$ by one and increasing the weight of the tree by $`w(x,u)w(x,v)w(u,v)`$.

Comment: In practice, we would wish to optimize the choice of neighbor of $`v`$ we choose to be adopted, as this changes the weight of the tree found.

### 2.1 The Adoption Network

###### Definitions 1

The deficit of vertex $`v`$ with respect to $`T`$ is $`deg_T(v)d(v)`$. Starting with a given tree, consider a sequence of adoptions $`\text{Adopt}(u_1,v_1),\text{Adopt}(u_2,v_2),\mathrm{}`$

* The sequence is legal if the precondition for each adoption is met.
* A sequence is feasible if, for each vertex, the decrease in its degree, i.e., its old degree minus its new degree, is at least its deficit.
* The cost of the sequence is $`_iw(u_i,v_i)`$.

The legal, feasible adoption sequences are precisely those that yield a tree meeting the degree constraints. The cost of a sequence is an upper bound on the resulting increase in the weight of the tree. Our goal is to find a feasible legal sequence of minimum cost. For brevity, we call such a sequence a minimum-cost sequence.

The problem reduces to a minimum-cost flow problem  in a flow network that we call the adoption network for $`T`$. The adoption network is defined as follows. Starting with $`G`$, replace each edge $`(u,v)`$ by two directed edges $`(u,v)`$ and $`(v,u)`$, each with cost $`w(u,v)`$ and infinite capacity. Assign each vertex a demand equal to its deficit.

A flow is an assignment of a real value (called the flow on the edge) to each edge of the network. The flow satisfies the following property, known as skew symmetry : $`f(u,v)=f(v,u)`$ for any edge $`(u,v)`$. For each vertex $`v`$, the surplus at $`v`$ is the net flow assigned to incoming edges minus the net flow assigned to outgoing edges. Only edges with positive flow are considered in computing the surplus flow into a vertex. A flow is legal if the surplus at each vertex is at most one less than its degree. A flow is feasible if the surplus at each vertex is at least its demand. The cost of the flow is the sum, over all edges, of the cost of the edge times the flow on the edge.

Since the demands are integers, there exists an integer-valued minimum-cost feasible flow . Assuming that each degree constraint is at least 1, there exists such a flow that is also legal. For brevity, we call such a flow a minimum-cost flow.

###### Lemma 1

The following statements are true:

1. The adoption sequences correspond to integer-valued flows. The correspondence preserves legality, feasibility, and cost.
2. The integer-valued flows correspond to adoption sequences. The correspondence preserves legality and feasibility; it does not increase cost.

Proof: Given a sequence of adoptions, the corresponding flow $`f`$ assigns a flow to each edge $`(u,v)`$ equal to the number of times $`u`$ adopts a neighbor of $`v`$. It can be verified that this correspondence preserves legality, feasibility, and cost.

Conversely, given an integer-valued flow $`f`$, modify it if necessary (by canceling flow around cycles) so that the set of edges with positive flow is acyclic. This does not increase the cost. Next, order the vertices so that, for each directed edge $`(u,v)`$ with positive flow, $`u`$ precedes $`v`$ in the order. Consider the vertices in reverse order. For each vertex $`u`$, for each edge $`(u,v)`$ with positive flow, have $`u`$ adopt $`f(u,v)`$ neighbors of $`v`$. Note that when an edge directed out of a vertex $`v`$ is processed by Adopt, $`v`$’s degree increases by one. Similarly when an edge directed into $`v`$ is processed, $`v`$’s degree decreases by one. The order imposed above processes all outgoing edges (in the flow) of a vertex before it processes any of its incoming edges. Hence during the course of the algorithm, $`v`$’s degree initially increases when its outgoing edges are processed and subsequently decreases as its incoming edges are processed. Therefore, when Adopt processes an edge $`(u,v)`$, the precondition imposed by the procedure that $`v`$’s degree be at least 2 always holds. It can be verified that the above sequence of adoptions preserves legality and feasibility and does not increase cost.

## 3 Polynomial-Time Algorithm

An acyclic, integer, minimum-cost flow can be found in polynomial time . The corresponding legal, feasible adoption sequence can be performed in polynomial time as described in the proof of the second part of Lemma 1. This gives a polynomial-time algorithm.

### 3.1 Optimality in Tree-Induced Metrics

The following lemma shows that this algorithm is optimal among algorithms that examine only the weights of edges of the given tree.

###### Lemma 2

Given a weighted graph $`G=(V,E)`$ and a spanning tree $`T`$ such that the weight of each edge in $`G`$ equals the weight of the corresponding path in $`T`$, a minimum-cost sequence of adoptions yields an optimal tree.

Proof: Fix an optimal tree. Note that the degree of $`v`$ in the optimal tree is at most $`d(v)`$; let it be $`\mathrm{d}^{}(v)`$. For each subset $`S`$ of vertices, let $`deg_T(S)`$ and $`\mathrm{d}^{}(S)`$ denote the sum of the degrees of vertices in $`S`$ in $`T`$ and in the optimal tree, respectively. Define a flow on the edges of $`T`$ as follows: for each edge $`(u,v)`$ in $`T`$, let $`f(u,v)=\mathrm{d}^{}(S_u)deg_T(S_u)`$, where $`S_u`$ is the set of vertices that are reachable from $`u`$ using edges in $`T`$ other than $`(u,v)`$. Note that $`f(u,v)=f(v,u)`$. Intuitively, a negative flow of $`x`$ units from $`u`$ to $`v`$ means that $`x`$ units of flow go from $`v`$ to $`u`$. This is known as the skew symmetry property of flows . Inductively it can be shown that for each vertex $`v`$, the net flow into it is $`deg_T(v)\mathrm{d}^{}(v)`$, so that the adoption sequence determined by the flow $`f`$ achieves a tree with the same degrees as the optimal tree.

We will show that the cost of the flow, and therefore the cost of the adoption sequence, is at most the difference in the weights of the two trees. This implies that the tree obtained by the adoption sequence is also an optimal tree.

To bound the cost of the flow, we claim that the flow is “necessary” in the following sense: for each edge $`(u,v)`$ in $`T`$, at least $`f(u,v)+1`$ edges in the optimal tree have one endpoint in $`S_u`$ and the other in $`VS_u`$. To prove this, let $`c`$ be the number of edges in the optimal tree crossing the cut $`(S_u,VS_u)`$. Note that $`deg_T(S_u)=2(|S_u|1)+1`$. Since the optimal tree is acyclic, the number of edges in the optimal tree with both endpoints in $`S_u`$ is at most $`|S_u|1`$. Thus $`\mathrm{d}^{}(S_u)2(|S_u|1)+c=deg_T(S_u)1+c`$. Rewriting gives $`c\mathrm{d}^{}(S_u)deg_T(S_u)+1=f(u,v)+1`$. This proves the claim.

To bound the cost of the flow, for each edge $`(u,v)`$, charge $`w(u,v)`$ units to each edge in the optimal tree crossing the cut $`(S_u,VS_u)`$. By the claim, at least the cost of the flow, plus the cost of $`T`$, is charged. However, since the cost of each edge in the optimal tree equals the weight of the corresponding path in $`T`$, each edge in the optimal tree is charged at most its weight. Thus, the total charge assigned to the edges is bounded by the weight of the optimal tree.

Note that given the exact degrees of the desired tree (for instance, if the degree constraints sum to $`2(|V|1)`$), the optimal flow in Lemma 2 can be computed in linear time.

### 3.2 Worst-Case Performance Guarantee

The next theorem establishes a worst-case performance guarantee for the algorithm in general graphs satisfying the triangle inequality.

###### Theorem 3

Given a graph $`G=(V,E)`$ with edge weights satisfying the triangle inequality, a spanning tree $`T`$, and, for each vertex $`v`$, a degree constraint $`d(v)2`$, the algorithm produces a tree whose weight is at most the weight of $`T`$ times

$$2\mathrm{min}\left\{\frac{d(v)2}{deg_T(v)2}:vV,deg_T(v)>2\right\}.$$

Proof: The increase in the cost of the tree is at most the cost of the best sequence. By Lemma 1, this is bounded by the cost of the minimum-cost flow. We exhibit a fractional feasible, legal flow whose cost is appropriately bounded. The minimum-cost flow is guaranteed to be at least as good.

Root the tree $`T`$ at an arbitrary vertex $`r`$. Push a uniform amount of flow along each edge towards the root as follows. Let $`p(v)`$ be the parent of each non-root vertex $`v`$. For a constant $`c`$ to be determined later, define

$$f(u,v)=\{\begin{array}{cc}c\hfill & \text{if }v=p(u)\hfill \\ 0\hfill & \text{otherwise.}\hfill \end{array}$$

The cost of the flow is $`c`$ times the weight of $`T`$. Let $`v`$ be any vertex. The surplus at $`v`$ is at least $`c(deg_T(v)2)`$. We choose $`c`$ just large enough so that the flow is feasible.

There are three cases. If $`deg_T(v)=1`$, the deficit at $`v`$ will be satisfied provided $`c1`$ and $`d(v)2`$. If $`deg_T(v)=2`$, the deficit at $`v`$ will be satisfied provided $`d(v)2`$. For $`deg_T(v)>2`$, the deficit will be satisfied provided

$$c\frac{deg_T(v)d(v)}{deg_T(v)2}=1\frac{d(v)2}{deg_T(v)2}.$$

Thus, taking

$$c=1\mathrm{min}\left\{\frac{d(v)2}{deg_T(v)2}:vV,deg_T(v)>2\right\}$$

gives the result.

## 4 Optimality of Performance Guarantee

In this section, we show that the worst-case performance guarantee established in Theorem 3 is the best obtainable.

###### Lemma 4

Consider an $`n`$-vertex weighted graph $`G`$ with a spanning tree $`T`$ such that the weight of each edge in $`T`$ is $`1`$ and the weight of each remaining edge is the weight of the corresponding path in $`T`$. If $`T`$ corresponds to a complete rooted $`(D1)`$-ary tree of depth $`k`$, then the weight of any spanning tree with maximum degree $`d`$ is at least the weight of $`T`$ times

$$2\frac{d2}{D2}o(1),$$

where $`o(1)`$ tends to $`0`$ as $`n`$ grows.

Proof: Fix any spanning tree $`T^{}`$ of maximum degree $`d`$. Let $`S_i`$ denote the vertices at distance at most $`i`$ from the root in $`T`$. The flow on the edges of $`T`$ corresponding to $`T^{}`$, as defined in Lemma 2, can be generalized to arbitrary cuts $`(VS,S)`$ in the tree, and it can be shown that the flow crossing this cut is at least $`deg_T(S)deg_T^{}(S)`$. For any $`i<k`$, the cut $`(VS_i,S_i)`$ is crossed by at least $`|S_i|(Dd)1`$ units of flow. Thus the total cost of the flow is at least $`_{i=0}^{k1}(|S_i|(Dd)1)`$. The cost of $`T`$ is $`|S_k|1`$, which can be written as $`_{i=0}^{k1}|S_{i+1}||S_i|`$. It can be verified that $`|S_{i+1}||S_i|=|S_i|(D2)+1`$. Hence the ratio of the cost of the flow to the cost of $`T`$ is at least

$$\frac{_{i=0}^{k1}(|S_i|(Dd)1)}{_{i=0}^{k1}(|S_i|(D2)+1)}.$$

On simplification, it can be seen that the ratio is at least $`(Dd)/(D2)o(1)`$. Since the ratio of the cost of $`T^{}`$ to the cost of $`T`$ is $`1`$ more than this, the result follows.

Next we observe that the $`d(v)2`$ constraint is necessary to obtain any constant performance guarantee:

###### Lemma 5

Consider an $`n`$-vertex weighted graph $`G`$ with a spanning tree $`T`$ such that the weight of each edge in $`T`$ is $`1`$ and the weight of each remaining edge is the weight of the corresponding path in $`T`$. If $`T`$ corresponds to a path of length $`n`$ with endpoint $`r`$, then the weight of any spanning tree in which each vertex other than $`r`$ has degree $`1`$ is at least the weight of $`T`$ times $`n/2`$.

The proof is straightforward.

## 5 Linear-Time Algorithms

Note that to obtain the worst-case performance guarantee a minimum-cost flow is not required. It suffices to find a feasible integer flow of cost bounded by the cost of the fractional flow $`f`$ defined in the proof of Theorem 3. We describe two methods to find such a flow, and to implement the corresponding sequence of adoptions, in linear time.

#### Algorithm 1:

The first algorithm exploits the special structure of the flow in the proof of Theorem 3 to construct an integral flow without solving the flow problem to optimality. Observe that the flow along each edge of the tree is the same ($`c`$ units). Since the graph satisfies the triangle inequality, flow along a path of more than a single edge can be replaced by a single edge that connects the end vertices of the path without increasing the total cost of the flow. This ensures that all of the fractional flows are sent from the sources to the destinations directly. The following greedy rounding scheme finds an integral flow that is no more expensive than the fractional flow.

Let $`f`$ be the fractional flow defined in Theorem 3. Modify $`f`$ by repeatedly performing the following short-cutting step: choose a maximal path in the set of edges with positive flow; replace the ($`c`$ units of) flow on the path by ($`c`$ units of) flow on the single new edge $`(u,v)`$, where the path goes from $`u`$ to $`v`$. Let $`q(u)`$ be the child of $`v`$ on the path. Stop when all paths have been replaced by new edges. This phase requires linear time, because each step requires time proportional to the number of edges short-cut.

In the resulting flow, the only edges with positive flow are edges from leaves of the (rooted) tree $`T`$ to interior vertices. Round the flow to an integer flow as follows. Consider each vertex $`v`$ with positive deficit, say $`D`$. Using a linear-time selection algorithm, among the edges $`(u,v)`$ sending flow to $`v`$, find the $`D`$ smallest-weighted edges. Assign one unit of flow to each of these $`D`$ edges. The resulting flow is integer-valued, feasible, legal, and has cost bounded by the cost of $`f`$. This phase requires linear time.

Assume that each vertex maintains a doubly linked list of its children. Given a pointer to any vertex, we can obtain its sibling in constant time. As adoptions are done, this list is maintained dynamically. Perform the adoptions corresponding to the flow in any order: for each edge $`(u,v)`$ with a unit of flow, have $`u`$ adopt the right sibling of $`q(u)`$ (in the original tree $`T`$). The tree remains connected because $`d(v)2`$, so at least one child of $`v`$ is not adopted.

#### Algorithm 2:

Consider the following restricted adoption network. Root the tree $`T`$ as in the proof of Theorem 3. Direct each edge $`(u,v)`$ of the tree towards the root. (Non-tree edges are not used.) Assign each edge a capacity of 1 and a cost equal to its weight. Assign each vertex a demand equal to its deficit.

We show below that an integer-valued minimum-cost flow in this network can be found in linear time. Because the fractional flow defined in the proof of Theorem 3 is a feasible legal flow in this network, the minimum-cost flow that we find is at least as good.

Find the flow via dynamic programming. For each vertex $`v`$, consider the subnetwork corresponding to the subtree rooted at $`v`$. Let $`C_j(v)`$ denote the minimum cost of a flow in this subnetwork such that the surplus at $`v`$ exceeds its demand $`D`$ by $`j`$, for $`j=0,1`$. Since the flow problem has been restricted to the tree, with a capacity constraint of 1 on all the edges, there is no need to consider flow surpluses greater than 1 at any vertex. For each child $`u`$ of $`v`$, let $`\delta (u)`$ denote $`w(u,v)+C_1(u)C_0(u)`$ — the additional cost incurred for $`v`$ to obtain a unit of flow along edge $`(u,v)`$. Let $`U_j`$ denote the $`D+j`$ children with smallest $`\delta (u)`$, for $`j=0,1`$. Then, for $`j=0,1`$,

$$C_j(v)=\underset{uU_j}{}\delta (u)+\underset{u}{}C_0(u).$$

Using this equation, compute the $`C_j`$’s bottom-up in linear time. The cost of the minimum-cost flow in the restricted network is given by $`C_0(r)`$, where $`r`$ is the root. The flow itself is easily recovered in linear time.

To finish, shortcut the flow as in the first phase of the previous algorithm and perform the adoptions as in the last phase of that algorithm.

## 6 Geometric Problems

Our general result has several implications for cases of particular distance functions where it is possible to give a priori bounds on the maximum degree of an MST. For the case of $`L_2`$ distances in the plane, there always is an MST of maximum degree 5 ; for the case of $`L_1`$ or $`L_{\mathrm{}}`$ distances there always exists a MST of maximum degree $`4`$ . Without using any specific structure of the involved distance functions, we note as a corollary:

###### Corollary 6

Let $`T_{\mathrm{min}}`$ be an MST and $`T_k`$ be a tree whose maximal degree is at most $`k`$. For $`L_1`$ or $`L_{\mathrm{}}`$ distances in $`\mathrm{IR}^2`$, we get a degree-3 tree $`T_3`$ with

* $`w(T_3)<\frac{3}{2}w(T_{\mathrm{min}})`$.

For the case of Euclidean distances in the plane, we get bounded degree trees that satisfy

* $`w(T_3)<\frac{5}{3}w(T_{\mathrm{min}})`$
* $`w(T_4)<\frac{4}{3}w(T_{\mathrm{min}})`$.

The latter two bounds are worse than those shown by Khuller, Raghavachari and Young  using the geometry of point arrangements. (It was shown that $`\frac{3}{2}`$ and $`\frac{5}{4}`$ are upper bounds.) We conjecture that the following are the optimal ratios:

###### Conjecture 7

For the case of Euclidean distances in the plane, we conjecture that there exist bounded degree trees that satisfy

* $`\frac{w(T_3)}{w(T_{\mathrm{min}})}\frac{\sqrt{2}+3}{4}1.103\mathrm{}`$
* $`\frac{w(T_4)}{w(T_{\mathrm{min}})}\frac{2\mathrm{sin}(\frac{\pi }{10})+4}{5}1.035\mathrm{}`$

For $`L_1`$ and $`L_{\mathrm{}}`$ distances in $`\mathrm{IR}^2`$, we conjecture

* $`\frac{w(T_3)}{w(T_{\mathrm{min}})}\frac{5}{4}`$

The worst examples we know (matching the ratios of Conjecture 7) are shown in Figure 2. (Note that the example for $`L_{\mathrm{}}`$ metric is obtained by rotating the arrangement in (c) by 45 degrees.)

### 6.1 Geometric Hamiltonian Paths

We conclude this paper by settling a question raised in , in the negative:

“In Euclidean graphs, perhaps a Traveling Salesman path of weight at most $`(2\epsilon )`$ times the minimum spanning-tree weight always exists and can be found in polynomial time.”

###### Theorem 8

For an arrangement of points in the plane with Euclidean distances, the ratio $`\frac{w(T_2)}{w(T_{\mathrm{min}})}`$ can be arbitrarily close to 2.

Proof: Let $`n`$ and $`k`$ be sufficiently large. Construct a point set as follows (see Figure 3):

Take base points at $`(0,0)`$ and $`(2n^k,0)`$.

For $`j=0,\mathrm{},k`$, add points as follows:

For $`i=1,\mathrm{},n^j`$

Add level $`j`$ point at $`((2i1)n^{kj},n^{kj})`$;

add base point at $`((2i1)n^{kj},0)`$.

The points at level $`j`$, i.e., at height $`n^{kj}`$, have nearest neighbors at distance at least $`n^{kj1}(n1)`$. To prove the lower bound, we draw a circle centered at each point at level $`j<k`$. For the points at level $`j`$, the radius of the circle is $`(n^{kj1}(n2))`$. The circles corresponding to two points do not intersect. Since each point has degree two in a Hamilton cycle, twice the sum of the radii of the circles gives us a lower bound on the length of the Hamilton cycle. This can be computed as follows (observe that we can always pick $`n2k`$).

$$2\underset{j=0}{\overset{k1}{}}n^j(n^{kj1}(n2))=2kn^{k1}(n2)2(k1)n^k.$$

Since no edge can have length more than $`2n^k`$, we conclude that no Hamilton path can have a weight smaller than $`2(k1)n^k2n^k=2(k2)n^k`$.

It can be verified that there is a tree of weight $`(k+3)n^k`$ that spans the points. Hence this is an upper bound on the weight of $`T_{\mathrm{min}}`$. It follows that $`\frac{w(T_2)}{w(T_{\mathrm{min}})}>\frac{2(k2)}{k+3}`$, which can be arbitrarily close to 2, concluding the proof.

The above class of examples establishes the same lower bound for $`L_1`$ and $`L_{\mathrm{}}`$ distances.

## Acknowledgements

We thank Joe Mitchell for establishing the transatlantic connection between the authors. We thank Chandra Chekuri for asking us about degree 3 trees in the $`L_1`$ metric.


text:
# The high-field phase diagram of the cuprates derived from the Nernst effect

## Abstract

Measurements of the Nernst signal in the vortex-liquid state of the cuprates to high fields (33 T) reveal that vorticity extends to very high fields even close to the zero-field critical temperature $`T_{c0}`$. In overdoped $`\mathrm{La}_{2\mathrm{x}}\mathrm{Sr}_\mathrm{x}\mathrm{CuO}_4`$ (LSCO) we show that the upper critical field $`H_{c2}(T)`$ curve does not end at $`T_{c0}`$, but at a much higher temperature. These results imply that $`T_{c0}`$ corresponds to a loss in phase rigidity rather than a vanishing of the pairing amplitude. An intermediate field $`H^{}(T)H_{c2}(T)`$ is shown to be the field scale for the flux-flow resistivity.

In the cuprate superconductors, strong fluctuations of the order parameter Fisher ; Blatter  and the extreme field scales make the task of constructing the experimental phase diagram a formidable challenge. While the important vortex-solid melting line $`H_m(T)`$ is now well-known Zeldov ; Schilling ; Liang , there is great uncertainty in the region above $`H_m`$. A key difficulty is the loss of long-range phase rigidity in the transition to the vortex-liquid state. In the vortex-liquid state, the flux-flow resistivity $`\rho `$ rises very rapidly to the (extrapolated) normal-state value Chien  even though substantial condensate strength and pairing amplitude remain. This makes $`\rho `$ an unreliable ‘diagnostic’ of the field-suppression of the pairing amplitude. By contrast, the Nernst and Ettinghausen effects can probe the existence of vorticity in the superfluid with undiminished sensitivity even in intense fields. From Nernst measurements in $`\mathrm{La}_{2\mathrm{x}}\mathrm{Sr}_\mathrm{x}\mathrm{CuO}_4`$ (LSCO) and $`\mathrm{YBa}_2\mathrm{Cu}_3\mathrm{O}_\mathrm{y}`$ (YBCO) in fields up to 33 T, we have uncovered an anomalous property of $`H_{c2}(T)`$ near their zero-field critical temperature $`T_{c0}`$. This anomaly is related to the existence of vortex-like excitations high above $`T_{c0}`$ Xu ; Wang , and is central to the key question of whether the Meissner transition in zero field is caused by the collapse of long-range phase coherence or the vanishing of the pairing amplitude (see also Corson et al.Corson ). In YBCO, previous Ettinghausen Palstra  and Nernst Ri ; Clayhold  measurements were performed in lower fields. Measurements on LSCO have been extended recently to high fields Capan .

In the Nernst experiment, vortices moving with velocity $`𝐯`$ down a thermal gradient $`T\widehat{𝐱}`$ generate a Josephson voltage which is observed as a transverse $`E`$-field $`E_y=Bv_x`$, with $`𝐁`$ the mean flux density. The vortex-Nernst effect is well-explored in low-$`T_c`$ superconductors. Figure 1b shows the Ettinghausen-Nernst effect in PbIn Vidal . The depinning of the vortex lattice by a supercurrent ($`𝐉\widehat{𝐱}`$) leads to a large transverse heat current $`J_y^h\widehat{𝐲}`$ that rises to a maximum and then decreases to zero at $`H_{c2}`$ (the weak fluctuation tail above $`H_{c2}`$ is not relevant here). By Onsager reciprocity, the Nernst signal $`e_yE_y/|T|`$ shares the same profile as $`J_y^h`$ in Fig. 1Onsager . To compare with theory, $`e_y`$ is divided by $`\rho (H)`$ to obtain the ‘transport line-entropy’ $`s_\varphi =\varphi _0e_y/\rho `$ ($`Ts_\varphi `$ is the heat energy carried by a unit-length vortex line). In PbIn, the linear decrease in $`s_\varphi `$ near $`H_{c2}`$ matches that of the magnetization $`M(T,H)=[H_{c2}(T)H]/1.16[2\kappa ^21]`$ ($`\kappa `$ is the Ginzburg-Landau parameter), as given by the Caroli-Maki expression Maki ; Vidal

$$s_\varphi (T,H)=\frac{\varphi _0}{T}|M(T,H)|L_D(T)$$
(1)
(here $`L_D(T)`$ decreases gradually with $`T`$ from 1 at $`T_{c0}`$).

We start with the overdoped regime, where the field profiles of $`e_y`$ most resemble those in PbIn. In Fig. 1a, we display curves of $`e_y`$ vs. field $`H`$ in Sample 1 \[$`\mathrm{La}_{2\mathrm{x}}\mathrm{Sr}_\mathrm{x}\mathrm{CuO}_4`$ (LSCO) with $`x`$ = 0.20 and $`T_{c0}`$ = 28 K\] (for experimental details see Ref. Wang ). At each temperature $`T`$, the signal rises steeply above $`H_m`$, attaining a prominent maximum before decreasing. The total data set defines experimentally the region in $`H`$ and $`T`$ where vorticity is strongly present. As in PbIn, the monotonic decrease at high fields reflects the field suppression of the condensate strength. At high fields, all the curves below 14 K are observed to follow a common curve towards zero (broken line). The trend implies that the low-$`T`$ traces are all wedged between the field axis and the broken line. Hence all the low-$`T`$ curves vanish at the intercept of the common curve with the field axis (45-50 T), which corresponds to $`H_{c2}(0)`$ .

Going to higher $`T`$, we immediately encounter an anomaly. Conventionally, the $`H_{c2}`$ line goes linearly to zero at $`T_{c0}`$. Hence $`e_y`$ ought to be finite in a field interval that $``$ 0 as $`TT_{c0}^{}`$. In sharp contrast, we find that, close to $`T_{c0}`$, the magnitude of $`e_y`$ remains large and nearly unchanged up to intense fields. There is no evidence for a field scale (flagged by $`e_y0`$) that decreases to zero regardless of how close we are to $`T_{c0}`$. To make this point quantitative, we convert $`e_y`$ to $`s_\varphi `$ using $`\rho (T,H)`$ measured in the same sample. The line entropy rises steeply to a sharp maximum before decreasing monotonically (Fig. 1c). The linear decrease at high fields, strikingly similar to that in PbIn, allows $`H_{c2}(T)`$ to be determined fit . These results show that $`H_{c2}(T)`$ decreases with rising $`T`$, but remains at extraordinarily high values as $`TT_{c0}`$ (see Fig. 2). \[We neglect a small correction from the negative Nernst coefficient $`\nu _N`$ -40 nV/KT of the holes in overdoped LSCO (Fig. 3 of Ref. Xu ). Accounting for the background increases $`s_\varphi `$ and $`H_{c2}(T)`$ by $`10\%`$. \]

A powerful way to summarize the anomalous trends described is to represent $`e_y(T,H)`$ as a contour plot in the $`T`$-$`H`$ plane (Fig. 2). In the vortex-solid phase below the melting line $`H_m`$ (lower white curve), $`e_y=0`$ because the vortices are pinned. As the vortex solid melts across $`H_m`$, $`e_y`$ increases steeply up to the ‘ridge’ which represents a field scale that we call $`H^{}(T)`$ (upper white curve). To the right of the ridge, $`e_y`$ falls monotonically. We focus on the two interesting regimes: the low-$`T`$ limit and the region near $`T_{c0}`$.

As $`T0`$, the two fields $`H_m(T)`$ and $`H^{}(T)`$ increase rapidly towards values approaching $`H_{c2}(0)`$ melting . In this large-$`H`$, low-$`T`$ region, the contours are nearly horizontal, as required by the asymptotic high-field behavior (broken curve in Fig. 1a). This implies that the $`H_{c2}`$ line is nearly horizontal below $``$12 K, as anticipated above. At temperatures $`TT_{c0}`$, the weak attenuation of $`e_y`$ to high fields corresponds to the nearly vertical orientation of the contours. Clearly, there is no evidence for a putative $`H_{c2}`$ line that terminates at $`T_{c0}`$. The values of $`H_{c2}(T)`$ derived from $`s_\varphi `$ are plotted in Fig. 2 as solid circles.

The anomalous features of the Nernst signal become more pronounced when we go to the underdoped regime. A major change from the overdoped case is apparent in the field profiles of $`e_y`$. We show in Fig. 3 results in underdoped YBCO (Sample 2, with $`y=6.50`$ and $`T_{c0}`$ = 50 K). As $`H`$ increases above $`H_m`$, $`e_y`$ rises rapidly, but attains a very broad maximum that extends undiminished to 30 T. In comparison with the curves in Fig. 1a, we see that $`e_y`$ in Sample 2 approaches zero only at fields considerably higher than 30 T. If we now convert $`e_y`$ to $`s_\varphi `$ using the measured $`\rho `$, the broad profile translates to a line entropy that extrapolates to zero at fields much higher than our maximum field. At 40 K, we may estimate the lower bound $`H_{c2}(T)`$60 T, which is very high for a temperature so close to $`T_{c0}`$ = 50 K. These features are common to all the underdoped cuprates we have investigated to date (for results in LSCO and Bi 2201, see Ref. Wang ).

The contour plot of $`e_y`$ for Sample 2 is displayed in the main panel Fig. 4. In comparison with Fig. 2, the melting line is now considerably lower relative to $`H^{}`$. Further, the contour lines in the region around $`T_{c0}`$ are more spread out. The Nernst signal retains considerable strength at 70 K, indicating that vorticity strongly persists up to 20 K above $`T_{c0}`$. The contour plots show rather clearly the continuity between the high-temperature fluctuation phase observed by Xu et al. Xu ; Wang  and the vortex-liquid state below $`T_{c0}`$. The 2 regimes smoothly merge together, and no phase boundary or ‘crossover’ line separating them is apparent in the $`T`$-$`H`$ plane. As mentioned above, the upper critical field is estimated to be above 60 T even at 40 K. The much higher field scale in underdoped YBCO is also apparent if we compare the contour plots in Figs. 2 and 4. In overdoped LSCO, the contours close at fields above 13 T (the field scale of the ‘island’) whereas, in underdoped YBCO, they do not close until $`H`$ exceeds 30 T.

For comparison, we also display the contours in near-optimal YBCO (Sample 3 with $`y`$ = 6.99 and $`T_{c0}`$ = 93 K). As shown by the vertical orientation of the contour lines from 85 to 95 K, as well as the field profiles in the insert in Fig. 3, $`e_y`$ is virtually field independent up to 14 T. Even for near-optimal YBCO, there is no evidence for an $`H_{c2}`$ line that decreases to zero at $`T_{c0}`$ (this calls into question the interpretation of ‘$`H_{c2}`$’ values inferred from bulk magnetization results Welp ).

The observation of vortex-like excitations high above $`T_{c0}`$ by Xu et al. Xu ; Wang  raised the question whether the cuprate superconducting transition at $`T_{c0}`$ actually corresponds to the loss of long-range phase rigidity, as opposed to the vanishing of the Gorkov pairing amplitude $`(𝐱,𝐱^{})`$ Emery . The present results provide strong support for the former. A prominent anomaly at temperatures near $`T_{c0}`$ is the remarkably weak field dependence of $`e_y`$ above the field scale $`H^{}(T)`$, which implies that the condensate strength extends undiminished to very high fields. In underdoped cuprates, this anomalous behavior persists to our lowest temperatures (Fig. 3). However, in overdoped samples, the profiles sharpen up below 25 K to reveal pronounced field suppression at high fields (Fig. 1a). Values of $`H_{c2}(T)`$ inferred from the line-entropy profiles show that the reason for the anomaly is that the $`H_{c2}`$ vs. $`T`$ line does not terminate at $`T_{c0}`$, but at a much higher temperature.

This seems to us to be compelling evidence that, at $`T_{c0}`$, the steep increase of $`\rho `$ and the collapse of the Meissner effect, correspond to a pre-emptive loss of phase rigidity (as has been advocated by Emery and Kivelson Emery  and others) rather than the vanishing of $`(𝐱,𝐱^{})`$. This scenario is qualitatively distinct from that in bulk low-$`T_c`$ superconductors. Although $`\rho `$ saturates rapidly to $`\rho _N`$ above $`T_{c0}`$ (in zero field), there remains sufficient condensate density to support vorticity. In Sample 1, the measured $`H_{c2}`$ curve implies that we need to apply fields above 30 T to suppress all vorticity near $`T_{c0}`$. The field scale is even higher in underdoped YBCO.

A hallmark of the vortex-liquid state in cuprates is that $`\rho `$ rapidly increases and ‘saturates’ to the extrapolated normal-state value ($`\rho _N`$) even when $`HH_{c2}(T)`$. To illustrate this, we plot the $`T`$ dependence of $`\rho `$ and $`e_y`$ in Sample 1 with $`H`$ fixed at 33 T in Fig. 1d (insert). Clearly, the large Nernst signal indicates that the region from 5 to 25 K is in the vortex-liquid state (see Fig. 2). However, between 15 and 30 K, $`\rho `$ virtually lies on the extrapolated $`\rho _N`$ curve. Hence reliance on $`\rho `$ alone would lead to a serious underestimate of ‘$`H_{c2}`$’.

In all cuprates studied, the ridge field $`H^{}(T)`$ line is ubiquitous. Over a broad range of $`T`$, the ratio $`\rho /\rho _N`$ has the value 0.3-0.4 on the curve $`H^{}`$ vs. $`T`$. $`H^{}(T)`$ serves as the boundary separating the low-field regime in which $`\rho `$ rapidly increases with $`H`$, from the high-field regime in which $`\rho `$ slowly asymptotes to the extrapolated $`\rho _N`$. Hence $`H^{}`$ represents an intrinsic field scale that controls dissipation in the vortex liquid state (attempts to find $`H_{c2}`$ using $`\rho `$ alone invariably turn up a curve akin to $`H^{}`$ rather than the higher $`H_{c2}`$). If we identify a length scale $`\xi ^{}(T)`$ by the ratio $`\xi ^{}(T)/\xi _0(T)\sqrt{H_{c2}(T)/H^{}(T)}`$, we find that at low $`T`$, $`\xi ^{}/\xi _01`$. As $`T`$ increases, however, the ratio rapidly increases (to 1.8 and 3.7 at 15 and 25 K, respectively). Since $`H^{}`$ determines the onset of strong flux-flow dissipation, $`\xi ^{}`$ represents a radius larger than the vortex core radius $`\xi _0`$ ($``$2.6 nm). Packing vortices closer than the outer radius $`\xi ^{}`$ leads to a crossover to the high-dissipation state in which $`\rho `$ asymptotes to $`\rho _N`$. However, as superfluidity and its associated 2$`\pi `$ phase winding survive right up to the inner radius $`\xi _0`$ Wen , the Nernst signal is observable to the higher field scale $`H_{c2}`$ ($`s_\varphi `$ is strongly attenuated in magnitude because of the reduced supercurrent for $`\xi ^{}<r<\xi _0`$fit . A length scale extending outside the vortex core has been observed in tunneling experiments Pan .

We thank S. Hannahs for valuable assistance at the National High Magnetic Field Lab. (a facility supported by NSF and the state of Florida). We acknowledge useful discussions with P.W. Anderson, A. J. Millis and Z. Y. Weng. This research is supported by a NSF grant (NSF-DMR 98-09483), the Office of Naval Research (N00014-01-0281), the New Energy and Industrial Tech. Develop. Org. (NEDO), Japan, the Natural Sciences and Engineering Research Council, and the Canadian Institute for Advanced Research.

Present address of ZAX: Department of Physics, Zhejiang University, Hangzhou, China


text:
# Azimuthal anisotropy of 𝐾_𝑆⁰ and Λ+Λ̄ production at mid-rapidity from Au+Au collisions at √𝑠__NN=130 GeV

## Abstract

We report STAR results on the azimuthal anisotropy parameter $`v_2`$ for strange particles $`K_S^0`$, $`\mathrm{\Lambda }`$ and $`\overline{\mathrm{\Lambda }}`$ at mid-rapidity in Au+Au collisions at $`\sqrt{s_{_{\mathrm{NN}}}}=130`$ GeV at RHIC. The value of $`v_2`$ as a function of transverse momentum of the produced particles $`p_t`$ and collision centrality is presented for both particles up to $`p_t3.0`$ GeV/c. A strong $`p_t`$ dependence in $`v_2`$ is observed up to 2.0 GeV/c. The $`v_2`$ measurement is compared with hydrodynamic model calculations. The physics implications of the $`p_t`$ integrated $`v_2`$ magnitude as a function of particle mass are also discussed.

Measurements of azimuthal anisotropies in the transverse momentum distribution of particles can probe early stages of ultra-relativistic heavy-ion collisions sorge99 ; sorge97 ; ollitrault92 . In high-energy nuclear collisions, the initial geometric anisotropy is established from the overlap between the colliding nuclei. The time necessary to build up this spatial anisotropy is believed to be short because the colliding nuclei are highly Lorentz contracted in the center-of-mass system and pass through each other at approximately the speed of light. During a $``$ 5–50 fm/c period, rescattering transfers the initial spatial anisotropy into a momentum anisotropy. This momentum anisotropy manifests itself most strongly in the azimuthal distribution of transverse momenta. The extent to which the initial spatial anisotropy is transformed to the measured momentum anisotropy depends on the initial conditions and the dynamical evolution of the system. In particular, anisotropy measurements for nucleus-nucleus collisions at RHIC energies may provide information about a partonic stage that may exist early in the collision evolution sorge99 ; gyulassy99 ; pasi01 ; gyulassy01 ; shuryak01 ; zlin01 .

The transverse momentum distribution of particles can be described in the form:

$$\frac{d^2N}{dp_t^2d\varphi }=\frac{dN}{2\pi dp_t^2}[1+2\underset{n}{}v_ncos(n\varphi )],$$
(1)
where $`p_t`$ is the transverse momentum of the particle, $`\varphi `$ is its azimuthal angle with respect to the reaction plane yingchao96 ; art98  and the harmonic coefficients, $`v_n`$, are anisotropy parameters. The second coefficient $`v_2`$ is called elliptic flow, where flow denotes collective behavior without necessarily implying a hydrodynamic limit. Recent experimental results from RHIC starflow1 ; starflow2 ; phenixv2qm01 ; phobosv2qm01  include measurements of $`v_2`$ as a function of collision centrality and $`p_t`$ for charged particles with $`p_t<2.0`$ GeV/c, and for identified charged pions, kaons and protons for $`p_t`$ up to $`0.8`$ GeV/c. The degree of the anisotropy transfer from position to momentum distribution depends on the density of the system during its evolution and the scattering cross sections of the particles involved (parton and/or hadron). As a result, recent theoretical work attempted to deduce the initial gluon density from partonic energy loss gyulassy01 , and the equation of state from hydrodynamic model calculations pasi01 ; shuryak01 .

Most of the anisotropic flow parameters measured to date are for non-strange particles starflow1 ; starflow2 ; ritter ; herrmann ; e87794 ; eos96 ; na4998 . Of the studies for identified strange particles starflow2 ; 4pi95 ; kaos98 ; eos98 ; e89501 ; e877la01 ; wa9897  most have been at much lower collision energies. Moreover, previous measurements of strange particle flow correspond to directed flow, i.e. the coefficient $`v_1`$. At the CERN SPS, quantitative differences between multi-strange baryons and non-strange hadrons were observed in transverse radial flow in Pb + Pb collisions at $`\sqrt{s_{_{\mathrm{NN}}}}`$ = 17 GeV na4497 ; wa9798 . A physical scenario in which multi-strange baryons do not participate in a common expansion and thus decouple early from the collision system due to their small hadronic cross sections, was proposed to explain this observation hecke98 . This explanation suggests that it may be possible to obtain insight into very early stages of the collisions by studying the elliptic flow of strange particles.

In this paper, we report the first measurement of the azimuthal anisotropy parameter $`v_2`$ for the strange particles $`K_S^0`$, $`\mathrm{\Lambda }`$ and $`\overline{\mathrm{\Lambda }}`$ from Au + Au collisions at $`\sqrt{s_{_{\mathrm{NN}}}}`$ = 130 GeV. Our measurement of $`v_2`$ for different centralities as a function of $`p_t`$ using the Solenoidal Tracker At RHIC (STAR) extends to a $`p_t`$ of about 3.0 GeV/c, much higher than previously measured for identified charged pions, kaons and protons starflow2 .

The STAR detector STAR , with its azimuthal symmetry and large acceptance, is ideally suited to measure elliptic flow. The detector consists of several sub-systems in a large solenoidal magnet. For collisions in its center, the Time Projection Chamber (TPC) measures charged tracks in the pseudo-rapidity range $`|\eta |<1.5`$ with $`2\pi `$ azimuthal coverage. During the year 2000 data taking the STAR magnet operated at a 0.25 Tesla field, allowing tracking of particles with $`p_t>0.075`$ GeV/c. A scintillator barrel surrounding the TPC, the Central Trigger Barrel (CTB), measures the charged particle multiplicity (for triggering) from within $`|\eta |<1`$. Two zero-degree calorimeters ZDC  located at $`\pm 18.25`$ m from the nominal interaction region, sub-tending an angle $`\theta <0.002`$ radians, primarily detect fragmentation neutrons. Two ZDCs in coincidence provide a minimum-bias trigger and the CTB is used for a central trigger. This analysis uses $`201\times 10^3`$ minimum-bias and $`180\times 10^3`$ central events.

We reconstructed both $`K_S^0\pi ^++\pi ^{}`$ and $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})p+\pi ^{}(\overline{p}+\pi ^+)`$ from their charged daughter tracks detected in the STAR TPC tpc . Using the energy loss of the charged tracks in the TPC gas, we select candidates for protons, anti-protons and pions. The mass and the kinematic properties of the neutral particle candidates are extracted from the decay vertex and daughter particle kinematics. Fig. 1 shows the invariant mass distributions for $`\pi ^+\pi ^{}`$ showing a $`K_S^0`$ mass peak and for $`p\pi ^{}`$ showing a $`\mathrm{\Lambda }`$ mass peak. The dashed lines are fits to the background and the signal. We determined that the background is dominated by combinatorial counts by rotating all positive tracks 180 degrees in the transverse plane and reconstructing the $`K_S^0`$ and $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})`$ decay vertices. This procedure destroys all real vertices in the TPC acceptance so that we can describe the combinatorial contribution to the invariant mass distributions. The observed masses, $`496\pm 8`$ MeV/$`c^2`$ for $`\pi ^+\pi ^{}`$ and $`1116\pm 4`$ MeV/$`c^2`$ for $`p\pi `$, are consistent with accepted values pdg  and the widths are determined by the momentum resolution of the detector. The particles used for the $`v_2`$ analysis are from the kinematic region of $`|y|1.0`$ and $`0.2p_t3.2`$ GeV/c for $`K_S^0`$ or $`0.3p_t3.2`$ GeV/c for $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$, where $`y`$ is the particle’s rapidity. No significant differences in elliptic flow are observed between $`\mathrm{\Lambda }`$ and $`\overline{\mathrm{\Lambda }}`$, so because of the limited statistics, $`\mathrm{\Lambda }`$ and $`\overline{\mathrm{\Lambda }}`$ are summed together.

We choose the requirements for $`K_S^0`$ and $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})`$ daughter candidates to maximize statistics. For $`K_S^0`$, we require the daughter tracks to have a distance-of-closest-approach (dca) to the collision vertex $`>`$ 1.0 cm. For the $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})`$ reconstruction, we choose pion-like tracks with a dca $`>`$ 1.5 cm and proton-like tracks with a dca $`>`$ 0.8 cm. We use the peak in the invariant mass distribution to measure the yield of $`K_S^0`$ or $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ particles for different values of $`\varphi `$ and $`p_t`$. Using the $`\varphi `$ bin center for the value of $`\varphi `$, we evaluate $`v_2`$ as a function of $`p_t`$ by calculating $`cos(2\varphi )`$ for different values of $`p_t`$. Using the yield to calculate $`v_2=cos(2\varphi )`$ enables us to measure elliptic flow for identified particles beyond the $`p_t`$ region where the identification of particles via their energy loss in the TPC gas fails starflow2 .

The real reaction plane is not known, but the event plane, an experimental estimator of the true reaction plane, can be calculated from the azimuthal distribution of tracks starflow1 . To determine the event plane, we select charged particle tracks with at least 15 measured space points, $`0.1<p_t2.0`$ GeV/c and $`|\eta |<1.0`$. We also require the ratio of the number of space points to the expected maximum number of space points for each track to be greater than 0.52, suppressing split tracks from being counted twice. Events are required to have a primary vertex within 75 cm longitudinally of the TPC center. These cuts are similar to those used in Ref. starflow1  and our analysis is not biased by them.

To avoid possible auto-correlations, tracks used for the $`K_S^0`$ or $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})`$ reconstruction are excluded from the set of tracks used to calculate the event plane. Typically this is done by measuring the azimuthal angle between a track and the event plane calculated from all other qualifying tracks within the same event. Then the contribution to $`v_2`$ from that track is calculated. In this analysis, where $`v_2`$ is not calculated on a particle by particle basis, all tracks that might be used for the reconstruction of $`K_S^0`$ or $`\mathrm{\Lambda }(\overline{\mathrm{\Lambda }})`$ are excluded from the event plane calculation. Only tracks with a dca $`<1.0`$ cm are used in the event plane calculation while the $`K_S^0`$ vertices don’t include these tracks. In the $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ analysis all proton-like tracks are excluded from the event plane calculation. A track is considered proton-like if its energy loss ($`dE/dx`$) is within three standard deviations of that expected for protons.

When the azimuthal anisotropy is evaluated via $`v_2=cos(2\varphi )`$, the observed $`v_2`$ must be corrected to account for the imperfect event plane resolution starno . This resolution is influenced by two factors that depend on centrality: the strength of the anisotropy signal and the number of tracks used for the event plane calculation. We estimate the resolution using the method of random subevents art98  and use the relative multiplicity, as in Ref. starflow1 , to measure the event centrality. The maximum resolution for the $`K_S^0`$ and $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ analysis is found to be $`0.681\pm 0.004`$ and $`0.582\pm 0.007`$ respectively and is reached in the centrality corresponding to 25–35% of the measured cross section. The relatively lower resolution for the $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ analysis is caused by the exclusion of proton-like tracks from the event plane calculation.

Elliptic flow as a function of transverse momentum for central and mid-central collisions calculated from $`201\times 10^3`$ minimum-bias and $`180\times 10^3`$ central events is shown in Fig. 2. The two particles show a similar $`p_t`$ dependence in the two centrality intervals. The $`p_t`$ dependence is stronger in more peripheral collisions than in the central collisions. A similar dependence was observed for charged particles in Au + Au collisions at the same RHIC energy starflow2 .

For this analysis, three main sources contribute to systematic errors in the measured anisotropy parameters: particle identification, background subtraction, and correlations unrelated to the reaction plane (non-flow) such as resonance decays, jets or Coulomb and Bose-Einstein correlations olliplb477 ; olliprc62 . The contribution from the first two sources is estimated by examining the variation in $`v_2`$ after changing several track and event cuts. We estimate that these effects contribute an error of less than $`\pm `$ 0.005 to $`v_2`$. The contribution to $`v_2`$ from non-flow effects, however, could be significant, especially in peripheral collisions. A previous study used the correlation of event plane angles from subevents to estimate the magnitude of these contributions starflowqm . Non-flow effects are assumed to contribute to the first and second harmonic correlations by similar amounts, so the magnitude of the first harmonic correlation sets a limit on the non-flow contributions to $`v_2`$. That study showed that the non-flow systematic errors for charged particles are typically +0 and -0.005, but are significantly larger in the more peripheral events where the error increases to +0 and -0.035 for the 58–85% most central events. These estimates are confirmed by measurements of $`v_2`$ using the 4th-order cumulant method, a method that is insensitive to non-flow effects but which leads to larger statistical errors aihong . We assume the systematic errors on $`v_2`$ for the neutral strange particles $`K_S^0`$ and $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ are similar to those found in the analysis of charged particles starflow2 .

To make a comparison with available hydrodynamic model calculations, we plot $`v_2(p_t)`$ for both $`K_S^0`$ and $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ from $`201\times 10^3`$ minimum-bias collisions in Fig. 3. The dashed lines represent the hydrodynamic model calculations pasi01  for (from top to bottom) pions, kaons, protons, and lambdas. Also shown in the figure is $`v_2(p_t)`$ for charged hadrons sqm2001proc . Within statistical uncertainty, the $`K_S^0`$ results are in agreement with the $`v_2`$ of charged kaons (not shown) starflow2 . We observe that $`v_2`$ for both strange particles increases as a function of $`p_t`$ up to about 1.5 GeV/c, similar to the hydrodynamic model prediction. In the higher $`p_t`$ region however ($`p_t2`$ GeV/c), the values of $`v_2`$ seem to be saturated. It has been suggested that the shape and height of $`v_2`$ above 2 – 3 GeV/c in a pQCD model is related to energy loss in an early, high-parton-density, stage of the evolution gyulassy01 .

The $`p_t`$ integrated anisotropy parameters for charged hadrons, $`K_S^0`$, and $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ from minimum-bias collisions are shown in Fig. 4. The integrated values of $`v_2`$ are calculated by parameterizing the yield with the inverse slope parameter of exponential fits to the $`K_S^0`$ or $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ transverse mass distributions sqm2001proc ; huithesis . The integrated $`v_2`$ is dominated by the region near the particle’s mean $`p_t`$ and is insensitive to the upper and lower bounds of the integration. Although the $`v_2(p_t)`$ of $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ is below the $`v_2(p_t)`$ of $`K_S^0`$ for most $`p_t`$, as shown in Fig. 3, the $`p_t`$ integrated $`v_2`$ values increase with the particle mass. This increase is partly due to the relatively higher mean $`p_t`$ of the $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$ compared to the $`K_S^0`$. In hydrodynamic models, although the spatial geometry of the pressure gradient and the resultant collective velocity are the same for all particles, massive particles tend to gain larger transverse momenta and so develop a larger elliptic flow. The hydrodynamic model calculations pasi01 , shown as a gray-band and central line, are, within errors, in agreement with this result. The width of the gray-band in Fig. 4 indicates the uncertainties of the model calculation, mostly due to the choice of the freeze-out conditions. The increase of $`v_2`$ with particle mass indicates that significant collective motion, perhaps established early in the collision, is an effective means to transfer geometrical anisotropy to momentum anisotropy. The nature of the particles during this process, however, whether parton or hadron, and the degree of thermalization for strange particles during the collective expansion remains an open issue.

In summary, we have reported the first measurement of the anisotropy parameter, $`v_2`$, for $`K_S^0`$ and $`\mathrm{\Lambda }+\overline{\mathrm{\Lambda }}`$, from Au + Au collisions at $`\sqrt{s_{_{\mathrm{NN}}}}=130`$ GeV. The $`v_2`$ values as a function of $`p_t`$ from mid-central collisions are higher at each $`p_t`$ than $`v_2`$ from central collisions. Hydrodynamic model calculations seem to adequately describe elliptic flow of the strange particles up to a $`p_t`$ of 2 GeV/c. For $`p_t`$ above 2 GeV/c, however, the observed $`v_2`$ seems to saturate whereas hydrodynamic models predict a continued increase with $`p_t`$. The $`p_t`$ integrated $`v_2`$ as a function of particle mass is consistent with a hydrodynamic picture where collective motion, established by a pressure gradient, transfers geometrical anisotropy to momentum anisotropy. Although the hadronic scattering cross sections of strange and non-strange particles may be different, we have yet to see deviations in the measured $`v_2`$ from hydrodynamic calculations at low $`p_t`$ for strange or non-strange particles. In a possible partonic phase prior to the hadronic epoch, the hadronic scattering cross sections for the final hadrons are not relevant. As such, if the elliptic flow of identified particles proves to be independent of their relative hadronic cross sections, it may be evidence that $`v_2`$ is established during a partonic phase.

Acknowledgments: We thank P. Huovinen for providing the results of the hydrodynamic model calculations. We wish to thank the RHIC Operations Group and the RHIC Computing Facility at Brookhaven National Laboratory, and the National Energy Research Scientific Computing Center at Lawrence Berkeley National Laboratory for their support. This work was supported by the Division of Nuclear Physics and the Division of High Energy Physics of the Office of Science of the U.S. Department of Energy, the United States National Science Foundation, the Bundesministerium fuer Bildung und Forschung of Germany, the Institut National de la Physique Nucleaire et de la Physique des Particules of France, the United Kingdom Engineering and Physical Sciences Research Council, Fundacao de Amparo a Pesquisa do Estado de Sao Paulo, Brazil, the Russian Ministry of Science and Technology and the Ministry of Education of China and the National Natural Science Foundation of China.


text:
# On the unitarity of S-matrix in 1-d case

## Abstract

It is shown that the scattering S-matrix is unitary even if the scattering potential $`U(x)`$ tends to different limits at $`x\pm \mathrm{}`$. This result is in contrast to the statements of some authors which argue that the unitarity may be broken in 1-d case if $`U(\mathrm{})U(+\mathrm{})`$. The mistake may result from a wrong normalization of the wave functions. This work may be considered as a comment to some of those works.

A fundamental property of the scattering S-matrix is its unitarity (c.f. e.g. Ref.). This is a consequence of the particle conservation law. However, there have appeared a number of works  where the S-matrix obtained or used is not unitary, and the authors have accepted this as a special feature of the problem owing to the different behavior of the potential at $`\pm \mathrm{}`$. This work may be considered as a comment on those works. It consists of two parts. First, we repeat a wrong derivation using a simple example and then prove the unitarity of the scattering matrix in the general case.

Let us consider a step potential of the form:

$`U(x)=\{\begin{array}{cc}0,x<0& \\ U_0,x>0& \end{array}`$ (3)
The general solution of the Schrödinger equation is well known:

$`\psi _l=Ae^{ikx}+Be^{ikx},x<0`$ (4)
$`\psi _r=Ce^{i\kappa x}+De^{\kappa x},x>0`$ (5)
where $`k=\sqrt{2E}`$ , $`\kappa =\sqrt{2(EU_0)}`$ and for simplicity we put $`\mathrm{}=m=1`$. Let us consider the most instructive case when $`E>U_0`$. The constants $`A,B,C,D`$ have to satisfy the relations $`\psi _l(0)=\psi _r(0)`$ and $`\psi _l^{}(0)=\psi _r^{}(0)`$. This yields

$$A+B=C+D,k(AB)=\kappa (CD)$$
(6)
The S-matrix components connect the incoming waves ($`Ae^{ikx}`$ and $`De^{\kappa x}`$) with the outgoing ones ($`Be^{ikx}`$ and $`Ce^{i\kappa x}`$). If one defines the S-matrix by the relation

$`\left(\begin{array}{cc}A& \\ D& \end{array}\right)=\stackrel{~}{S}\left(\begin{array}{cc}C& \\ B& \end{array}\right)`$ (11)
it will be of the form

$$\stackrel{~}{S}_{11}=\frac{2\kappa }{k+\kappa },\stackrel{~}{S}_{22}=\frac{2k}{k+\kappa },\stackrel{~}{S}_{12}=\stackrel{~}{S}_{21}=\frac{k\kappa }{k+\kappa }$$
(12)
The quantity $`R=\stackrel{~}{S}_{12}^2=\stackrel{~}{S}_{21}^2`$ coincides with the reflection coefficient, but $`\stackrel{~}{S}_{11}`$ and $`\stackrel{~}{S}_{22}`$ have no physical meaning. They are not transmission amplitudes, at least because $`\stackrel{~}{S}_{11}^{}{}_{}{}^{2}+\stackrel{~}{S}_{12}^{}{}_{}{}^{2}1`$. Therefore, $`\stackrel{~}{S}`$ is not the scattering matrix of our problem.

Let us now derive the correct form of the S-matrix. We shall do this for an arbitrary potential of the form:

$$U(x\mathrm{})=0,U(x\mathrm{})=U_0$$
(13)
Therefore, the asymptotes at $`\pm \mathrm{}`$ are of the form (4). However, for calculating the scattering S-matrix components these functions have to be normalized to unit flux, not to unit density! Only in this case the S-matrix components correspond to reflection and transmission amplitudes (the corresponding coefficients are ratios of particle fluxes). So, we write, instead of (4), the following functions:

$$\psi _L=A_1\frac{e^{ikx}}{\sqrt{k}}+B_1\frac{e^{ikx}}{\sqrt{k}},x\mathrm{}$$
(14)
$$\psi _R=A_2\frac{e^{i\kappa x}}{\sqrt{\kappa }}+B_2\frac{e^{\kappa x}}{\sqrt{\kappa }},x\mathrm{}$$
(15)
In (7) and (8) $`k`$ and $`\kappa `$ coincide with the particle velocities at infinity because we put $`\mathrm{}=m=1`$. The functions $`\psi _L`$ and $`\psi _R`$ are asymptotes of one and the same exact solution. Hence, the following relations between coefficients $`A_i`$ and $`B_i`$ exist :

$$A_2=\alpha A_1+\beta B_1,B_2=\beta ^{}A_1+\alpha ^{}B_1$$
(16)
Since the wave functions are now normalized to unit fluxes, the particle conservation law (or the continuity condition for the current) reads

$$|A_1|^2|B_1|^2=|A_2|^2|B_2|^2$$
(17)
Substituting here $`A_2`$ and $`B_2`$ from (9) leads to

$$|\alpha |^2|\beta |^2=1$$
(18)
One has from (9):

$$A_1=\frac{1}{\alpha }A_2\frac{\beta }{\alpha }B_1,B_2=\frac{\beta ^{}}{\alpha }A_2+\frac{1}{\alpha }B_1$$
(19)
This means, that the scattering matrix has the form:

$`\widehat{S}=\left(\begin{array}{cc}1/\alpha & \beta /\alpha \\ & \\ \beta ^{}/\alpha & 1/\alpha \end{array}\right)`$ (23)
It is easy to check using (11), that this matrix is unitary. The reflection and transmission coefficients are respectively

$$R=|\beta /\alpha |^2=|\beta /\alpha ^{}|^2,T=|\alpha |^2$$

and $`T+R=1`$ as it should be.

It is instructive to note in conclusion, that the problem of normalization does not appear in the typical 3-d situations where the scattering potential vanishes at infinity. Nevertheless it is helpful sometimes to use the normalization to unit flux as e.g. in the Born approximation (c.f. Ref.  § 126).


text:
# References

Specificity and stability in topology of protein networks

Sergei Maslov1 and Kim Sneppen2

1 Department of Physics, Brookhaven National Laboratory,
Upton, New York 11973, USA
2 Department of Physics, Norwegian University of Science and Technology, N-7491 Trondheim, Norway.

Abstract:

Molecular networks guide the biochemistry of a living cell on multiple levels: its metabolic and signalling pathways are shaped by the network of interacting proteins, whose production, in turn, is controlled by the genetic regulatory network. To address topological properties of these two networks we quantify correlations between connectivities of interacting nodes and compare them to a null model of a network, in which al links were randomly rewired. We find that for both interaction and regulatory networks, links between highly connected proteins are systematically suppressed, while those between a highly-connected and low-connected pairs of proteins are favored. This effect decreases the likelihood of cross talk between different functional modules of the cell, and increases the overall robustness of a network by localizing effects of deleterious perturbations.

With the growth of experimental information about basic biochemical mechanisms of life, molecular networks operating in living cells are becoming better defined. Direct physical interactions between pairs of proteins form one such network. It serves as a backbone for functional and structural relationships among its nodes and defines pathways for the propagation of various signals such as phosphorylation and allosteric regulation of proteins. The information about specific binding of proteins to each other has recently grown by an unprecedented amount as a result of high throughput two-hybrid experiments . The production and degradation of proteins participating in the interaction network is controlled by the genetic regulatory network of the cell formed by all pairs of proteins in which the first protein directly regulates the abundance of the second. The majority of known cases of such regulation happens at the level of transcription, in which a transcription factor positively or negatively regulates the RNA transcription of the controlled protein. The large scale structure of both these networks is characterized by a high degree of interconnectedness, where most pairs of nodes are linked to each other by at least one path. One may wonder how such a heavily intertwined and mutually dependent dynamical system can perform multiple functional tasks, and remain stable against deleterious perturbations.

We analyzed the topological properties of interaction and transcription regulatory networks in yeast Saccharomyces cerevisiae, which at present is perhaps the best characterized model organism. The interaction network used in this work consists of 4549 physical interactions between 3278 yeast proteins as measured in the most comprehensive two-hybrid screen of yeast proteins , while the genetic regulatory network is formed by 1289 directed positive or negative direct transcriptional regulations within a set of 682 proteins as listed in the YPD database . The protein interaction network is a representative of the broad class of scale-free networks  in which the number of nodes with a given number of neighbors (connectivity) $`K`$ scales as a power law $`1/K^\gamma `$. In our case the histogram of connectivities can be fitted by a power law with $`\gamma =2.5\pm 0.3`$ for $`K`$ ranging from 2 to about 100 . A small part of the protein interaction network, formed by proteins known to be localized in the nucleus and to interact with at least one other nuclear protein, was visualized (Fig. 1). One striking feature of this graph is the abundance of highly connected proteins that are mostly connected to those with low connectivity, and thus well separated from each other.

To test for correlations in connectivities of nodes for each of the above two networks we calculated the likelihood $`P(K_0,K_1)`$ that two proteins with connectivities $`K_0`$ and $`K_1`$ are connected to each other by a link and compared it to the same quantity $`P_r(K_0,K_1)`$ measured in a randomized version of the same network. In this “null model” network all proteins have exactly the same connectivity as in the original one, while the choice of their interaction partners is totally random. The transcription regulatory network is naturally directed, while the network of physical interactions among proteins in principle lacks directionality. However, for poorly understood reasons the two-hybrid experimental data have a significant asymmetry between baits and preys, with bait hybrids being more likely to be highly connected than their prey counterparts. This can be seen e.g. in the fact that average connectivity of baits with at least one interaction partner is close to 3, whereas the same quantity measured for preys is only 1.8. Since each reported interaction involves one bait and one prey protein, this asymmetry needs to be taken into account when constructing an uncorrelated “null” model for the interaction network. For this purpose in our randomization procedure we would treat the two-hybrid data as a directed network with an arrow on each edge pointing out from bait to prey hybrid. Randomized versions of these two networks were constructed by randomly reshuffling links, while keeping the in- and out-degree of each node constant. A convenient numerical algorithm performing such randomization consists of first randomly selecting a pair of directed edges A$``$B and C$``$D The two edges are then rewired in such a way that A becomes connected to D, while C to B. However, in case if one or both of these new links already exist in the network this step is aborted and a new pair of edges is selected. This last restriction prevents the appearance of multiple edges connecting the same pair of nodes. A repeated application of the above rewiring step leads to a randomized version of the original network. Multiple sampling of randomized networks allowed us to calculate both the average expectation and the standard deviation for any particular property of the random network.

Correlations in connectivities manifest themselves as systematic deviations of the ratio $`P(K_0,K_1)/P_r(K_0,K_1)`$ from 1. We calculated this ratio for interaction (Fig. 2A) and regulatory (Fig. 2B) networks, with $`K_0`$ and $`K_1`$ being the total number of interaction partners of two interacting proteins (for the interaction network), and out- and in-degrees of two nodes connected by a directed edge 0$``$1 (for the regulatory network). Thus by the very construction $`P(K_0,K_1)`$ is symmetric for the physical interaction network but not for the regulatory network. We also estimated the statistical significance $`Z(K_0,K_1)`$ of the above deviations in the interaction (Fig. 2C) and regulatory (Fig. 2D) networks, by dividing each observed deviation from the null model by the standard deviation in multiple realizations of a randomized network. The combination of these two plots reveals the regions on the $`K_0K_1`$ plane, where connections between proteins in the real network are significantly enhanced or suppressed, compared to the null model. In particular red regions in the upper left and the lower right corners reflect the tendency of highly connected nodes (hubs) to associate with nodes of low connectivity, while the blue/green region in the upper right corner reflects the reduced likelihood that two hub centers are directly linked to each other. One should also note a prominent feature on the diagonal of the Fig. 2A and 2C corresponding to an enhanced affinity of proteins with between 4 and 9 interaction partners to physically interact with each other. This feature can be tentatively attributed to the tendency of members of multi-protein complexes to interact with other proteins from the same complex. The above range of connectivities thus correspond to a typical number of direct interaction partners of a protein in a complex. When we checked for interactions between proteins in this range of connectivities we found 39 pairs of interacting proteins to belong to the same complex in a recent high throughput study , which is 4 times more than one would expect to find by pure chance alone.

To further quantify and compare correlation patterns in interaction and regulatory networks we calculated the average connectivity $`K_1`$ of nearest neighbors of a node, as a function of its own connectivity $`K_0`$ (Fig. 3A). In order to simplify the comparison between two networks here we characterize each node in the regulatory network by its total number of neighbors $`K=K_{in}+K_{out}`$. For both interaction and regulatory networks the average connectivity $`K_1`$ shows a gradual decline with $`K_0`$, which can be fitted with a power law $`K_11/K_0^{0.6\pm 0.1}`$ over approximately two decades. This observation gives an additional credit to the affinity between correlation patterns in these two protein networks visible in Fig. 2. It was recently found  that the internet, defined as the set of interconnected routers, in addition to a scale-free distribution of node connectivities similar to the protein interaction network, is characterized by the same correlation pattern between connectivities of neighboring nodes: $`K_11/K_0^{0.5}`$. This extends by one step an intriguing similarity in the topology of these networks of completely different nature.

For the scale-free physical interaction network we also plotted the probability distribution of the nearest neighbor connectivity $`K_1`$, measured separately for nodes with small connectivity $`K_03`$, and for those with large connectivity $`K_0100`$ (Fig. 3B). In the absence of correlations this conditional probability does not depend on $`K_0`$, and is proportional to $`K_1/K_1^\gamma 1/K_1^{1.5}`$. This uncorrelated form holds approximately true for neighbors of a protein with low connectivity. It is only violated at the far tail of the distribution due to an excess likelihood of it being connected to a protein with very high connectivity, as was mentioned above. On the other hand, the distribution of connectivities $`K_1`$ of neighbors of highly connected proteins scales as $`1/K_1^{2.5}`$ and thus differs from that of lowly connected ones by a factor of $`1/K_1`$.

When analyzing molecular networks one should consider possible sources of errors in the underlying data. Two-hybrid experiments give rise to false positives of two kinds. In one case the interaction between proteins is real but it never happens in the course of the normal life cycle of the cell due to spatial or temporal separation of participating proteins. In another case an indirect physical interaction is mediated by one or more unknown proteins localized in the yeast nucleus. Reversely, in a high throughput two-hybrid screens one should expect a sizeable number of false negatives. Primarily a binding may not be observed if the conformation of the bait or prey heterodimer blocks relevant interaction sites or if the corresponding heterodimer altogether fails to fold properly. Secondly, 391 proteins out of the potential 5671 baits in  were not tested as possible bait hybrids because they were found to activate transcription of the reporter gene in the absence of any prey proteins.

Unlike for the interaction network, our data for the genetic regulatory network do not come from a single large scale project. Instead, they represent a collection of numerous experiments performed by different experimental techniques in different labs. Therefore, it is not feasible even to list possible sources of errors present in such a diverse data set. In particular one should worry about a hidden anthropomorphic factor present in such a network: some proteins just constitute more attractive subjects of research and are, therefore, relatively better studied than others. One should also note that the transcription regulation network is only a subset of a larger genetic regulatory network, which in addition to transcriptional regulation includes translational regulation, RNA editing, etc. An encouraging sign was that when we separately analyzed the set representing the current knowledge  about this later more complete network, consisting of 1750 genetic regulations among 848 proteins we reproduced all of our empirical results for the transcriptional network.

The observed suppression of connections between nearest neighbors of highly connected proteins is consistent with compartmentalization and modularity characteristic of control of many cellular processes . In fact, it suggests the picture of functional modules of the cell organized around individual hubs. To further test the extent of modularity of hubs and their immediate neighborhood in each network we selected 15 highest connected nodes. To provide an unbiased sample of hubs from the point of view of in and out connectivity half of those nodes were selected as the highest out-degree hubs (8 baits with $`K_{bait}90`$ for the interaction network and 7 nodes with $`K_{out}34`$ for the regulatory network), while half were the highest in-degree hubs (7 preys with $`K_{prey}20`$ for the interaction network and 8 nodes with $`K_{out}8`$ for the regulatory network). In agreement with the correlation properties described above, direct connections between hubs were significantly suppressed. In the interaction network we observed 20 links between different hubs in this group, which is significantly below $`56\pm 7.5`$ links in the randomized network. In the transcription regulatory network there were 16 links between hubs in real network as opposed to $`35\pm 6.5`$ in its randomized version. Not only direct links between hubs are suppressed in both studied networks, but hubs also tend to share fewer of their neighbors with other hubs, thereby extending their isolation to the level of next-nearest neighbor connections. The total number of paths of length 2 between the set of 15 hubs in the interaction network is equal to 418, whereas in the null model we measured this number to be $`653\pm 56`$. Similarly, for the transcriptional network the number of paths of length 2 is equal to 186 in the real network, whereas from the null model one expects it to be $`262\pm 30`$. Since the number of paths of length 2 between a pair of proteins is equal to the number of their common interaction partners one concludes that both the hub node itself and its immediate surroundings tend to separate from other hubs, reinforcing the picture of functional modules clustered around individual hubs.

A further implication of the observed correlation is in the suppression of the propagation of deleterious perturbations over the network. It is reasonable to assume that certain perturbations such as e.g. significant changes in the concentration of a given protein (including its vanishing altogether in a null-mutant cell) with a ceratin probability can affect its first, second, and sometimes even more distant neighbors in the corresponding network. While the number of immediate neighbors of a node is by definition equal to its own connectivity $`K_0`$, the average number of its second neighbors, given by $`K_0(K_11)_{K_0}`$, is sensitive to correlation patterns of the network. Since highly connected nodes serve as powerful amplifiers for the propagation of deleterious perturbations it is especially important to suppress this propagation beyond their immediate neighbors. It was argued that scale-free networks in general are very vulnerable to attacks aimed at highly connected nodes . The anticorrelation presented above implies a reduced branching ratio around these nodes and thus provides a certain degree of protection against such attacks. This may be the reason why the correlation between the connectivity of a given protein and the lethality of the mutant cell lacking this protein is not particularly strong .

It is feasible that molecular networks in a living cell have organized themselves in an interaction pattern that is both robust and specific. Topologically the specificity of different functional modules can be enhanced by limiting interactions between hubs and suppressing the average connectivity of their neighbors. We have seen that such correlation pattern appears in a similar way in two different layers of molecular networks in yeast, and thus presumably is a universal feature of all molecular networks operating in living cells.