When RD-WZQ model is used, the quantization MSE for U_m is given by 14

d m RD - WZQ = λ m 2 - 2 r m ,

where λ_m = var(U_m| Y) and var(·|·) denotes the conditional variance. The optimal solution to the WZ-TC problem under RD-WZQ model is given by the following theorem.

Theorem 1 Given jointly Gaussian X and Y as defined above, and a total bit budget of B bits, if each transform coefficient U_m, m = 1, . . ., M ₁, where U = T ^T X, is quantized by an RD optimal WZ quantizer which uses Y as decoder side-information, then the transform T which minimizes E X - X ^ 2 is the CKLT of X given Y, and the number of bits allocated to quantizing U_m is

ρ m ( λ , B , N ) = 1 2 log 2 λ m d * ( λ , B , N ) m = 1 , … , N 0 m = N + 1 , … , M

where N ≤ M is the largest integer for which λ_m ≥ d^∗ ( λ , B, N), m = 1, . . ., N. The quantization MSE of the mth coefficient is

d m = λ m 2 - 2 ρ m = d * ( λ , B , N ) , m = 1 , … , N λ m , m = N + 1 , … , M

and the overall MSE is

E X - X ^ 2 = N d * ( λ , B , N ) + ∑ m = N + 1 M 1 λ m .

Proof 1 Directly follows from [ 8 , Section III-B].

Note that RD-WZQ model implies infinite-dimensional VQ of each coefficient and hence the above MSE is the OPTA in the Gaussian WZ-TC problem.

The asymptotically (in rate) optimal solution to the WZ-TC problem under SWC-HRSQ model is given by the following theorem.

Theorem 2 Let X, Y, and B be as in Theorem 1. If each transform coefficient U_m , m = 1, . . ., M ₁ , where U = T ^T X, is quantized by a high-rate scalar quantizer and the quantizer output is encoded by a SW code which uses Y as decoder side-information, then the transform T which asymptotically minimizes E X - X ^ 2 is the CKLT of X given Y and the bit allocation is given by (2). Furthermore, the asymptotically optimal quantizer for U_m, m = 1, . . ., N is a uniform quantizer with step-size Δ = ( 2 π e ) d * ( λ , B , N ) . The resulting quantization MSE is

E X - X ^ 2 = π e 6 N d * ( λ , B , N ) + ∑ m = N + 1 M 1 λ m .

Proof 2 See Section "Proof of Theorem 2" in Appendix 1.

The WZ quantizers with vector-valued decoder side-information as considered in Theorem 1 are difficult to design in practice. However, the following theorem establishes that when CKLT is used and RD-WZQ model applies for quantization of coefficients, a linear transformation of the side-information vector can be used to convert the vector side-information problem into an equivalent scalar side-information problem. Furthermore, [ 5 , Section 6] shows that this result applies in asymptotic sense to the SWC-HRSQ model as well.

Theorem 3 Let the mean-zero vectors X ∈ ℝ M 1 and Y ∈ ℝ M 2 be jointly Gaussian, and let T be the CKLT of X given Y. Suppose that transform coefficients U_m , m = 1, . . ., M ₁, where U = T ^T X, are each compressed by an RD optimal WZ quantizer relative to decoder side-information Y. Then, the minimum MSE (MMSE) estimate ũ_m( y ) = E{U_m| y} of U_m given Y = y is a sufficient statistic for decoder side-information for quantizing U_m .

Proof 3 See Section "Proof of Theorem 3" in Appendix 1.

Wyner-Ziv transform coding is a special case of more general MT-TC where two or more terminals apply TC to their respective inputs and transmit the quantized outputs to a single decoder which exploits the inter-source correlation to jointly reconstruct all the sources. In this case, the problem is to optimally allocate a given bit budget among all the terminals such that the total MSE is minimized. However, the closed-from solution to this problem appears difficult, due to the inter-dependence of the encoders in different terminals. An iterative descent algorithm is given in 8 for solving the Gaussian MT-TC problem. Given a total bit-budget, the bit-rate of the system is incremented by a small amount in each iteration, and the optimal WZ-TC for each terminal is determined by fixing the encoders of all other terminals and considering their outputs as decoder-side information. The solution that gives the MMSE is accepted and the iterations are repeated until the total bit-budget is exhausted. While this algorithm, referred to as the DKLT algorithm, is guaranteed to converge to at least a locally optimal solution, there is no tractable way to implement the quantizers implied by the final solution since it is not practical to optimize a set of near-optimal WZ quantizers in each iteration of this algorithm. Note also that, DKLT requires joint decoding of two vector sources.

We note that the optimal solution to the unconstrained problem defined in (8) corresponds to the MMSE solution of the constrained problem over the set of rate-tuples S = B 1 ′ , B 1 ″ , B 2 : B 1 ′ ∈ 0 , B , B 1 ″ ∈ 0 , B , B 2 ∈ 0 , B , B 1 ′ + B 1 ″ + B 2 ≤ B . This set is shown in Figure 2. One approach to locating the MMSE solution is to search over an appropriately discretized grid of points inside . As we will see, even though an exhaustive search on a fine grid can be prohibitively complex, a much simpler constrained tree-search algorithm exists which can be used to locate the required solution with a very high probability.

The proposed algorithm is a generalization of a class of bit-allocation algorithms in which a small fraction ΔB of the total bit-budget B is allocated to the "most deserving" quantizer among a set of quantizers in an incremental fashion, until the entire bit-budget is exhausted [ 4 , Section 8.4]. Unfortunately, this type of a greedy search cannot guarantee that the final solution is overall optimal and can yield poor results in our problem where the bit allocation among three sets of dependent quantizers must be achieved. On the other hand, if the increment ΔB is chosen small enough, a near-optimal solution can be found by resorting to a tree-search. Even though a full tree-search is intractable, a simple algorithm referred to as the (M, L)-algorithm 18 exists for detecting the minimum cost path in the tree with a high probability. We use this insight to formulate a tree-search algorithm for solving the unconstrained bit allocation problem, in which a set of constrained bit allocation problems are solved in each iteration.

In order to describe the proposed tree-search algorithm in detail, let ΔB be the incremental amount of bits to be allocated in each step of the search, where 0 < ΔB ≪ B. The algorithm is initialized by setting B 1 ′ , B 1 ″ , B 2 = 0 , 0 , 0 , i.e., the origin in Figure 2. Now if we are to allocate ΔB bits to only one of the three transform codes T 1 ′ , T 1 ″ or T ₂, then there are three possible choices for the rate-tuple B 1 ′ , B 1 ″ , B 2 , namely (ΔB, 0, 0), (0, ΔB, 0), and (0, 0, ΔB). For each of these choices, we can explicitly solve the constrained bit allocation problem as described in the previous section and find the MMSE solution. Each of these candidate solutions can be viewed as a node in a tree as shown in Figure 3. The root node of the tree corresponds to a SP-DTC of rate 0, and a node in the first level of nodes obtained in the first iteration of the algorithm corresponds to a SP-DTC of rate of ΔB bits per source pair (X ₁, X ₂). In the second iteration of the algorithm, we allocate ΔB more bits to each of the three candidate SP-DTCs (but to one of the SP-DTCs at a time) in the first level of nodes. Note that, for each SP-DTC we can allocate ΔB bits in three different ways, i.e., ΔB bits can be added to either B 1 ′ , B 1 ″ , or B ₂. This requires the solution of three constrained bit allocations problems for each of the 3 nodes in the first-level. As a result, the tree will be extended to a second level of 3² nodes, in which each node corresponds to a SP-DTC of 2ΔB bits per source-pair, as shown in Figure 3. We can repeat this procedure, allocating ΔB bits to each of terminal node of the tree in a given iteration, until all B bits are exhausted. After the final iteration, the tree would consist of L = ⌈B/ΔB⌉ levels with 3 ^L nodes in the last level (terminal nodes). Each terminal node corresponds to a candidate SP-DTC of rate B, and rate-tuples of these SP-DTCs lie on the plane B 1 ′ + B 1 ″ + B 2 = B in Figure 2. The MMSE terminal node of the tree is the optimal solution to the unconstrained bit allocation problem, provided that the latter solution is on the search-grid. If ΔB is chosen small enough, then we can ensure that the optimal solution is nearly on the search-grid. Suppose that, in each iteration, the algorithm saves the MSE of the solution to the constrained bit-allocation problem associated with each node. In theory, the optimal solution can be found by an exhaustive tree-search, using the MSE of a node [given by (22)] as the path-cost. In order to practically implement the tree-search, we use the (M, L)-algorithm, in which the parameter M can be chosen to reduce the complexity at the expense of decreased accuracy (i.e., the probability of detecting the lowest cost path in the tree). In the (M, L) algorithm [ 18 , p. 216] for a tree of depth L, one only retains the M best (lowest MSE) nodes in each iteration. When M = 1 we have a completely greedy search. On the other hand when M_n = 3 ⁿ in the nth iteration, we have a full-tree search which has a complexity that grows exponentially with the iteration number. When M (1 ≤ M ≤ 3 ^L ) is a prescribed constant, the complexity is linear in M, independent of the number of iterations. In obtaining the simulation results presented in Section 4, M = 27 and L = 135 (ΔB = 0.2) were found to be sufficient to obtain near optimal results. For example, it was observed that even for M = 81 and L = 405, nearly the same result was obtained.

We compute the rate-pairs (R ₁, R ₂) achievable with a SP-DKLT code for a given a total MSE D, by fixing R ₁ (or R ₂) and then searching for minimum R ₂ (or R ₁) required to achieve the MSE D [given by (6)]. The rate-pairs achievable for D = 0.01 with SP-DKLT coding of vectors from source model A with ρ = 0.9, γ = 0.9, are plotted in Figure 4. These values of ρ and γ result in a source cross-covariance matrix with the largest element 0.9. In Figure 4, the curve "SP-DKLT (X ₁ split)" corresponds to a system in which the input to terminal 1 (which applies source splitting) is X ₁ as shown in Figure 1. The curve "SP-DKLT (X ₂ split)" corresponds to a system in which the input to the terminal 1 is X ₂. Note that the two curves are not symmetric in rates and they coincide if R ₁ and R ₂ are inter-changed in one of the curves. This is because X ₁ and X ₂ have different auto-covariance matrices, and hence inter-changing the rates is equivalent to inter-changing the terminals. Importantly, this result indicates that when the two sources are not statistically identical, which source is chosen for splitting does not affect the SP-DKLT performance. Table 1 lists the best bit allocations found by the tree-search algorithm for SP-DKLT codes shown in Figure 4. Note that, for the same (B ₁, B ₂), the rate-split between B 1 ′ and B 2 ″ when X ₁ is applied to the terminal 1 is not identical to that when X ₂ is applied to the terminal 1.

Figure 4 also shows the rate region achievable if each source is independently compressed using the KLT (i.e., only intra-vector correlation is utilized), labeled IKLT, and the OPTA lower bound for distributed TC predicted by the iterative DKLT algorithm 8 . The performance of both distributed and non-distributed TC of source model A degrades as ρ decreases, since both auto- and cross-covariance matrices of X ₁ and X ₂ are functions of ρ. Next consider the source model B for which the lowest achievable rate-pairs corresponding to D = 0.005 are plotted in Figure 5. The source parameter α = 0.32 results in auto-covariance matrices whose largest off-diagonal element is 0.9. Also recall that θ is the largest element in the source cross-covariance matrix. Note that changing θ only affects the cross-covariance matrix, and hence has no effect on the best achievable rates for independent coding of the two sources. On the other hand, as θ increases, the rates achievable with distributed coding do improve. Since in source model B, the two sources are statistically similar, the curves in Figure 5 are symmetric in rates and the optimal bit allocation does not depend on which source is chosen for splitting.

The RD performance in Figures 4 and 5 indicate that SP-DKLT codes can significantly outperform IKLT codes at all rates when there is sufficient correlation between the two distributed sources. The performance of SP-DKLT coding necessarily approaches OPTA (DKLT) bound when either R ₁ or R ₂ is sufficiently high. That is, the terminal with the higher bit rate can independently transform code its input with negligible distortion, and the other terminal can then apply WZ-TC at the minimum rate achievable with "almost unquantized" decoder side-information. However, for both source-models the rate-region achievable with source-splitting is strictly inside that of DKLT. In other words, there are some rate-pairs inside the DKLT rate-region for a given MSE D, which cannot to be achieved by a SP-DKLT code. A closely related issue is that, for a range of values of (R ₁, R ₂), the sumrate R ₁ + R ₂ of SP-DKLT codes remains constant and reaches its minimum. For example, it can be seen from Figure 4 and Table 1 that the sum-rate is about 4.125 bits when the rate of X ₁ is in the range 1.375 - 2.5 bits/sample. From Table 1, it can be seen that when the sum-rate is greater than its minimum value, the optimal SP-DKLT code approaches a WZ transform code, i.e., no source-splitting occurs. This situation, which also exists in Figure 5, suggests that optimal SP-DKLT codes for the sum-rate at which a given D can be achieved, are equivalent to time-sharing 12 of two "corner points". Figure 6 illustrates this situation for optimal SP-DKLT codes at D = 0.005 for source model B (θ = 0.9 in Figure 5). It should however be noted that, unlike source-splitting, time-sharing between the two terminals requires synchronization of their encoders 11 .

In this section, we focus on the practical design of SP-DKLT codes for a given pair of rates (R ₁, R ₂) based on both scalar and block-quantization. RD-WZQ model used in Section 3.1.1 implies infinite block-length WZ-VQ of each coefficient. A practically realizable approach to block WZ quantization is SWC-TCQ 20 . Experimental results obtained with LDPC codes of block length up to 10⁶ bits and TCQs up to 8,192 states are presented in 20 for quadratic Gaussian WZ quantization, which indicate that performance very close to the theoretical limit can be achieved with SWC-TCQ. Motivated by these results, we aim to implement SP-DTCs which can approach theoretical performance predicted in Section 3.1.1 using TCQ and SW coding for encoding transform coefficients. However, the SWC-TCQ design procedure followed in 20 is to first design a TCQ whose MSE satisfies a constraint (by choosing a sufficiently high rate) and then to estimate the output conditional entropy (which is the target rate of the SW code) of the resulting TCQ. This is sufficient for verifying the achievable rate pairs for a given MSE which is the goal of 20 . Our problem is different in that the rate of the SW code is specified by the solution to the bit-allocation problem and our goal is to design a TCQ which minimizes the MSE, subject to a constraint on the output conditional entropy. This requires an alternative formulation of the design procedure, which we refer to as CEC-TCQ. In previous work on non-distributed quantization, entropy constrained TCQ (EC-TCQ) has been investigated in 21 22 23 24 . CEC-TCQ is a modification of EC-TCQ in 21 22 to accommodate block SW-coding of the TCQ output relative to a decoder side-information sequence. Our formulation of CEC-TCQ follows the supersetentropy formulation of EC-TCQ in 22 .

Suppose that a sequence of source samples {U_n ∈ℝ} has to be quantized, given that the sequence {Y_n ∈ℝ}, is available at the decoder as side-information, where n = 1, 2, . . . denotes the discrete-time. Similar to an ordinary TCQ 25 , a CEC-TCQ uses a size 2 R TCQ + 1 scalar codebook to quantize the input sequence U ₁, U ₂, . . ., into a R _TCQ bits/sample output sequence Û ₁, Û ₂, . . . . However, CEC-TCQ output satisfies the additional property that the conditional entropy H (Û_n|Y_n ) = E{- log₂ P (Û_n|Y_n )} ≤ R for some given R. It follows from 9 that if the CEC-TCQ output is SW-coded relative to the decoder side-information sequence {Y_n } then LR bits are sufficient to (almost) losslessly transmit a sequence of L source samples as L → ∞. The optimal CEC-TCQ minimizes the MSE E{(U_n - Û_n )²}, subject to the constraint E{- log₂ P (Û_n|Y_n )} ≤ R, or equivalently, minimizes the Lagrangian

J = lim L → ∞ ∑ n = 1 L E U n - Û n 2 + β E - log 2 P Û n | Y n

where β > 0 is the Lagrange multiplier. This implies that, given a specific sequence of input samples u ₁, u ₂, . . ., the CEC-TCQ encoder should use the Viterbi algorithm based on the path-cost function

∑ n u n - û n 2 + β E - log 2 P û n | Y n .

In a rate R _TCQ TCQ, each codeword c_k , k = 1 , … , 2 R TCQ + 1 , in the codebook is labeled with R _TCQ-bit binary string 25 . Let b_i , i = 1 , … , 2 R TCQ be these binary labels. Then, to compute (25), the cost βE{- log₂ P (b_i|Y)} must also be stored for each binary label, where Y is the random variable representing the decoder side-information. For a fixed value of β, we can use a slight modification of the algorithm in 21 for optimizing the CEC-TCQ codebook, by replacing the codeword entropies E{- log₂ P(b_i )} by the conditional entropies E{- log₂ P (b_i|Y)}, and by using a training sequence of (U_n, Y_n ) pairs. In order to approximate the expectations by sample averages computed from training data, the side-information variable Y is discretized to Ŷ ∈ η 1 , … η Y , where is a large enough positive integer. Then E{- log₂ P (B = b_i |Y)} can be approximated by

Ĥ b i = - ∑ k = 1 Y log 2 P B = b i | Ŷ = η k P Ŷ = η k | B = b i

where B ∈ b 1 , … , b 2 R TCQ is the binary-labeled output of the TCQ. Given a TCQ code-book, the probabilities P(B = b_i |Ŷ = η_k ) and P(Ŷ = η_k |B = b_i ) can be estimated using the training set [it is sufficient to estimate p_{i, k} = P(B = b_i, Ŷ = η_k ), k = 1 , … Y , i = 1 , … , N ]. To complete the design, it is necessary to search for the value of β for which E{- log₂ P (Û_n |Y_n )} ≈ R by repeating the codebook optimization for an appropriately chosen sequence of β values.

For block WZ-code designs, the transforms and the bit allocations are found by RD-WZQ model (Section 3.1.1) and WZ quantizers are implemented using CEC-TCQ followed by binary SW coding. More specifically, the rate found by the bit allocation algorithm for each transform coefficient is used as the conditional entropy constraint in the design of a CEC-TCQ for that coefficient. As described in Section 2.3, the CEC-TCQ designs are based on scalar-side information obtained by a linear transform of the vector side-information at the decoder, see Theorem 3. All CEC-TCQ designs are based on the 8-state trellis used in JPEG2000 [ 26 , Figure 3.16]. For trellis encoding and decoding, a sequence length of 256 source samples has been used. For design and testing quantizers, sample sequences of length 5 × 10⁵ have been used. Since, the main focus this paper is the design of transforms and the quantizers, we assume ideal SW coding of the binary output of each CEC-TCQ, so that our results do not depend on any particular SW coding method. In a practical implementation (e.g., 20 ), near optimal performance can be obtained by employing a sufficiently long SW code (note that sequence length for SW-coding can be chosen arbitrary larger than the sequence length used for TCQ encoding). This type of coding is well suited for applications such as distributed image compression, where the coding is inherently block-based.

We also consider SP-DKLT code designs based on scalar quantization. In this case, the transforms and bit allocations are found by using the SWC-HRSQ model (Section 3.1.2). While it is possible to use the step-size predicted by SWC-HRSQ model to design uniform quantizers, we found that such quantizers in reality do not satisfy the required entropy constraint at lower rates. We instead use conditional entropy constrained scalar quantizers (CEC-SQ), designed by modifying the algorithm in 27 to accommodate a conditional entropy constraint similar to CEC-TCQ approach above.

The reconstruction signal-to-noise ratio (RSNR) [with MSE as given by (6)] of SP-DKLT code designs for source model B is shown in the rows labeled Design in Table 2 where SP-DKLT/CECSQ and SP-DKLT/CECTCQ refer to scalar quantization and TCQ based designs respectively. The rows labeled Analytical show the performance predicted by the SWC-HRSQ and RD-WZQ models upon which the transforms and bit-allocations are based (note however that the performance predicted by SWC-HRSQ model is not necessarily an upper-bound for CEC-SQ designs which are not constrained to be uniform quantizers). We compare the performance of our SP-DKLT code designs with IKLT codes for which the bit-allocations are obtained by using either entropy coded high-rate quantization model (for scalar quantizer design) [ 4 , Section 9.9] or RD-optimal quantization model (for block quantizer design) [ 16 , Section 10.3.3] for Gaussian variables. The IKLT codes with scalar quantization have been implemented by using entropy constrained scalar quantizers (EC-SQ) while those with block quantizers have been implemented by using EC-TCQ 23 , where we assume ideal entropy coding of the quantizer outputs. In Table 2, IKLT/ECSQ and IKLT/ECTCQ respectively refer to these designs.

From a practical view point, the use of DCT instead of KLT is interesting 26 . We therefore consider the design of SP-DTCs based on the DCT, referred to as source-split distributed DCT (SP-DDCT) codes. Since DCT is a fixed transform, we only need to optimize the bit allocations. To do this, we assume that DCT is approximately a decorrelating transform for Gaussian vectors 4 . Then, the bit-allocations given by (9), (14), and (20) are still valid provided that the eigenvalues in these expressions are replaced by the variances of the corresponding DCT coefficients. The rest of the design procedure is the same as that with SP-DKLT. The RSNR of DCT based designs are presented in Table 3 (again, the performance predicted by SWC-HRSQ model is not an upper bound for corresponding practical codes). The results show that, for this particular source model, even the scalar quantization based SP-DDCT codes outperform TCQ-based IDCT codes.

The optimality of CKLT is proved in 5 . To prove the optimality of the bit allocation, consider high-rate scalar quantization of the transform coefficient U_m and ideal SW coding of the quantizer output Û_m at the rate r_m = H (Û_m | Y) bits/sample, m = 1, . . ., M ₁, where H (·|·) denotes the conditional entropy 16 . In this case, the asymptotically optimal scalar quantizer for each coefficient is known to be uniform 5 . For high-rate uniform quantization, H (Û_m| Y) ≈ h(U_m| Y) - log₂(Δ _m ), where Δ _m is the quantizer step-size and h(·|·) is the conditional differential entropy [ 16 , Section 8.3]. Since the conditional variance of the transform coefficient U_m , given the side-information Y is E U m 2 | Y = λ m , and (U_m , Y) are jointly Gaussian, it follows that h(U_m| Y) = (1/2) log₂ (2πeλ_m ) 16 , and hence Δ m = 2 π e λ m 2 - r m . Therefore, the MSE of high-rate uniform quantization followed by ideal SW coding of U_m is

d m SWC - HRSQ = Δ m 2 12 = π e 6 λ m 2 - 2 r m .

Since (27) and (1) are the same within a constant factor of πe/6 (which is identical for all transform coefficients), it is easy to verify that the optimal bit-allocation solution under SWC-HRSQ model is also given by (2). However, the MSE of mth coefficient in this case is

d m = π e 6 λ m 2 - 2 ρ m = π e 6 d * λ , B , N , m = 1 , … , N λ m , m = N + 1 , … , M

and hence the overall MSE is given by (5). Furthermore,

Δ m = Δ = 2 π e d * λ , B , N

for m = 1, . . ., N.

For jointly Gaussian and mean-zero X and Y, there exists a matrix A and a mean-zero Gaussian vector W ₁ independent of Y such that X = AY + W ₁, where A = ∑ X Y ∑ Y - 1 , ∑ W 1 = ∑ X | Y = T Λ T T , and Λ is the diagonal matrix of eigenvalues of ∑ _X|Y . Furthermore U = T ^T AY + W ₂, where W ₂ = T ^T W ₁ is an uncorrelated Gaussian vector since ∑ W 2 = T T ∑ W 1 T = Λ (note that for CKLT, T ^-1 = T). Therefore ∑ _U|Y = Λ. The MMSE estimate of U given Y is Ũ = E{U|Y} = T ^T E{X | Y} = T ^T AY. Thus, U = Ũ + W ₂ and ∑ _U|Ũ = Λ, and it follows that U_i is independent of Ũ_j if j ≠ i and var(U_i| Ũ_i ) = var(U_i | Y), where var(·|·) denotes the conditional variance. Now, since h(U_i| Y) = h(U_i| Ũ_i ), we conclude that Ũ_i is a sufficient statistic 16 for decoder side-information Y in WZ quantization of U_i .