Hervé J Égou

Christine Guillemot

Robust multiplexed codes for compression of heterogeneous data

Keywords: source coding, variable length codes, data compression, data communication, entropy codes

Compression systems of real signals (images, video, audio) generate sources of information with different levels of priority which are then encoded with variable length codes (VLC). This paper addresses the issue of robust transmission of such VLC encoded heterogeneous sources over error-prone channels. VLCs are very sensitive to channel noise: when some bits are altered, synchronization losses can occur at the receiver. This paper describes a new family of codes, called multiplexed codes, that confine the de-synchronization phenomenon to low priority data while asymptotically reaching the entropy bound for both (low and high priority) sources. The idea consists in creating fixed length codes for high priority information and in using the inherent redundancy to describe low priority data, hence the name multiplexed codes. Theoretical and simulation results reveal a very high error resilience at almost no cost in compression efficiency.

I. INTRODUCTION

Entropy coding, producing variable length codewords (VLC), is a core component of any data compression scheme. Unlike fixed length codes (FLCs), variable length codes are designed to exploit the inherent redundancy of a symbol distribution. The main drawback of VLCs is their high sensitivity to channel noise: when some bits are altered by the channel, synchronization losses can occur at the receiver, the position of symbol boundaries are not properly estimated, leading to dramatic symbol error rates. Soft VLC decoding ideas, exploiting residual source redundancy (the so-called "excess-rate"), have also been shown to reduce the "de-synchronization" effect as well as the residual symbol error rates. These ideas rely on capitalizing on source coder suboptimality, by exploiting inner codeword redundancy and exploiting correlation within the sequence of symbols (inter symbol dependency). Models incorporating both VLC-encoded sources and channel codes (CC) have also been considered [START_REF] Murad | Joint source-channel decoding of variable length encoded sources[END_REF], [START_REF] Demir | Joint source-channel coding for variable length codes[END_REF], [START_REF] Bauer | Turbo fec/vlc decoding and its application to text compression[END_REF]. The authors in [START_REF] Murad | Joint source-channel decoding of variable length encoded sources[END_REF] derive a global stochastic automaton model of the transmitted bitstream by computing the product of the separate models (Markov source (MS), source coder (SC) and channel coder (CC)). The authors in [START_REF] Guyader | Joint sourcechannel turbo decoding of entropy coded sources[END_REF] push further the above idea by designing an iterative estimation technique alternating the use of the three models (MS, SC, and CC). The above methods mostly consist in re-augmenting the redundancy of the bitstream, by introducing an error correcting code or dedicated patterns in the chain. Also, the decoding algorithms, often relying on MAP estimators, have a rather high complexity.

The synchronization recovery capability (or selfsynchronization property) of a VLC, which represents the error-resiliency of the code, has also been considered as a performance and design criterion of the code in addition to its compression performance. The characterization of the synchronization capability of VLC has thus been extensively studied in the research community [START_REF] Wei | On the characterization of statistically synchronizable codes[END_REF], [START_REF] Capocelli | On the characterization of statistically synchronizable variable length codes[END_REF]. This property is often characterized by an error span [START_REF] Maxted | Error recovery for variables length codes[END_REF], also called the mean error propagation length (MEPL) [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF]. The authors in [START_REF] Maxted | Error recovery for variables length codes[END_REF] have developed a state model for synchronization recovery which is suitable for analyzing the performance of various codes with respect to error recovery. Effort has then also naturally been dedicated to the design of self-synchronizing codes. In [START_REF] Ferguson | Self-synchronizing huffman codes[END_REF], the authors show that some Huffman codes, tailored to a given source, may contain a codeword that can serve as a "re-synchronization" point. However, the existence of this codeword depends on the source statistics. Statistically synchronizable codes containing a synchronization sequence are also described in [START_REF] Wei | On the characterization of statistically synchronizable codes[END_REF] and [START_REF] Capocelli | On the construction of statistically synchronizable variable length codes[END_REF]. Self-synchronizing Huffman codes are also proposed in [START_REF] Lam | Self-synchronizing variable-length codes for image transmission[END_REF]. The design is governed by a trade-off between redundancy and synchronization recovery capability. E.g., the authors in [START_REF] Takashima | Error states and synchronization recovery for variable length codes[END_REF] and [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF] search to construct codes with minimum redundancy which allow to have short mean error propagation length (MEPL). In the same spirit, reversible VLCs (RVLC) [START_REF] Takishima | Reversible variable length codes[END_REF], [START_REF] Wen | Reversible variable length codes for efficient and robust image and video coding[END_REF], [START_REF] Bauer | Turbo fec/vlc decoding and its application to text compression[END_REF] have been designed specifically to fight against desynchronizations. Variable-tofixed length codes [START_REF] Tunstall | Synthesis of noiseless compression codes[END_REF] [START_REF] Savari | Generalized tunstall codes for sources with memory[END_REF] mainly designed for compression purposes are also well-suited for robust transmission. Both self-synchronizing codes and variable-to-fixed codes only allow for synchronization in terms Levenshtein distance sense [START_REF] Levenshtein | Binary codes capable of correcting deletions, insertions and reversals[END_REF]. However, strict-sense synchronization [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF], i.e. synchronization in the Hamming distance sense, is required by applications, such as image and video compression.

Here, we consider the design of a new family of codes, called "multiplexed codes", that allow to control (even avoid) the "de-synchronization" phenomenon of VLC encoded sources, while still allowing to reach asymptotically the entropy bound and to rely on simple decoding techniques. The principle underlying these codes builds upon the idea that compression systems of real signals generate sources of information with different levels of priority (e.g. texture and motion information for a video signal). This idea underlies unequal error protection (UEP) techniques, which allocate different levels of protection to the different types of information, either to signal frequency bands [START_REF] Tanabe | Subband image coding using entropycoded quantization over noisy channels[END_REF], to bit planes [START_REF] Ruf | Rate-distortion performance for joint source and channel coding of images[END_REF], or to quantization layers [START_REF] Cheung | Joint source-channel coding of scalable video over noisy channels[END_REF]. In the wake of this idea, we consider two sources, a high priority source and a low priority source referred to as S H and S L in the sequel. The codes designed are such that the risk of "de-synchronization" is confined to the low priority information. The idea consists in creating a FLC for the S H source, and in exploiting the inherent redundancy to represent or store information of the low priority source S L . Hence, the source S H inherits some of the properties of FLCs such as random access in the data stream and strictsense synchronization. It is shown that these codes allow to almost reach the entropy bound.

The rest of the paper is organized as follows. Section II introduces the notations we use in the sequel. Section III outlines the principle of multiplexed codes, discusses their compression properties and describes the first coding algorithm in the case of two sources. The approach relies on a mapping of the low priority flow of data into a sequence of k-valued variables, where k depends on the realization of the high priority sequence of symbols. A first mapping, optimal in the sense that it allows to reach the entropy bound is presented. This mapping is based on an Euclidean decomposition of a long integer, which can be computed with a hierarchical (dyadic) approach with a close to linear computational complexity. The computational cost can be further reduced (but this time at the expense of a small excess-rate), by performing the decomposition on a constrained set of k-valued variables. This amounts to consider a constrained partition of the set of fixed length codewords into equivalence classes expressed in terms of prime-valued variables. The resulting loss in compression depends on the distance between the respective probability density functions (pdf) of the constrained set of k-valued variable and of the source S H . A heuristic method for selecting the set of k-valued variables (i.e., the partition) that would best approximate the pdf of the source S H is described in section VI. A second class of multiplexed codes for which the cardinalities of the different equivalence classes are constrained to be integer powers of two is also presented. Such codes can be built in a very straightforward manner from any VLC code. The compression performance is then the same as the one achieved by the VLC considered. The impact of channel noise, considering a binary symmetric channel, on the overall source distortion is analyzed in section VII. The choice of the parameters of the algorithm are discussed in section VIII and simulation results are provided in section IX.

II. PROBLEM STATEMENT AND NOTATIONS

Let S H = (S 1 , . . . S t , . . . S KH ) be a sequence of source symbols of high priority taking their values in a finite alphabet A composed of |A| symbols, A = {a 1 , . . . a i , . . .}. Note that, in the following, we reserve capital letters to random variables, and small letters to values of these variables. Bold face characters will be used to denote vectors or sequences. The stationary probability function density of the source S H is denoted µ = (µ 1 , . . . µ i , . . . µ |A| ), where µ i stands for the probability that a symbol of S H equals a i . Let h be the stationary entropy of the source per symbol, given by h = -|A| i=1 µ i log 2 (µ i ). Let S L be a sequence of source symbols of lower priority taking their values in a finite alphabet. We assume that the realization s L of this source has been preencoded into a bitstream b = (b 1 , . . . b r , . . . b KB ) with a VLC coder (e.g. Huffman or arithmetic coder). The problem addressed here is the design of a family of joint codes for the two sources S H and S L that would guarantee the strict-sense synchronization of the high priority source S H at almost no cost in terms of compression efficiency.

class C i codeword x symbol a i |C i | probability µ i index q i C 1 000 a 1 3

III. MULTIPLEXED CODES

A. Principle and definitions

Let us denote X the set of binary codewords of length c. This set of |X | = 2 c codewords is partitioned into |A| subsets, denoted C i for i = 1 . . . |A|, called equivalence classes associated to symbols a i of the alphabet A, as shown in Table I. Note that codebook partitioning is sometimes referred to as binning. The condition c ≥ log 2 (|A|) is required to have at least 1 codeword per subset. Each equivalence class C i is a set of |C i | codewords. In the following, the integer |C i | denotes the cardinality of the equivalence class and is such that

|A| i=1 |C i | = |X |.
A symbol S t = a i of the flow S H can be encoded with any c-bit codeword x belonging to the equivalence class C i (see example of Table I). Hence, each codeword S t can be mapped into a pair (C i , Q i ) of two variables denoting respectively the equivalence class and the index of the codeword in the equivalence class C i . The variable Q i is an |C i |-valued variable, taking its value between 0 and |C i | -1 (see Table I) and representing the inherent redundancy of the c-bits FLC. This redundancy will be exploited to describe the lower priority flow of data. Let us denote n t the cardinality |C i | of the equivalence class associated to the realization S t = a i . Then, to the realization s H = s 1 . . . s t . . . s KH of the sequence of high priority symbols one can associate a sequence q 1 . . . q t . . . q KH of n t -valued variables that will be used to describe jointly the high and low priority data flows (i.e., s H and s L ).

Definition 1: a multiplexed code is defined as a set of pairs (a i , q i ), where a i is a symbol value belonging to the alphabet A (on which the high priority source is quantized) and where q i is the value of a variable Q i denoting the index of the codeword in the equivalence class.

The corresponding encoder is the function which maps every couple (a i , q i ) into an c-bits fixed length codeword x as ( a i , q i ) ⇆ x.

(

) 1
Example 1: let A = {a 1 , a 2 , a 3 , a 4 , a 5 } be the alphabet of the source S H with the stationary probabilities given by µ 1 = 0.4, µ 2 = 0.2, µ 3 = 0.2, µ 4 = 0.1 and µ 5 = 0.1. Note that this source has been considered in [START_REF] Maxted | Error recovery for variables length codes[END_REF] and [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF]. Table I gives an example of partitioning of the set X into the 5 equivalence classes associated to the alphabet symbols. This partitioning reveals 5 |C i |-valued variables. Note that, in Table I, the codes are ranked in the lexicographic order. Using this order, each equivalence class can be fully defined by, 1) the first codeword within this equivalence class, 2) the cardinality of this equivalence class. Thus, the lexicographic order permits a compact representation of the multiplexed code table. However, note that the method is not restricted to this order.

B. Conversion of the lower priority bitstream

It appears from above that the design of multiplexed codes relies on the choice of the partition C = {C 1 , . . . C |A| }. The encoding then proceeds with the conversion of the lower priority bitstream into a flow of variables taking their values in the different sets of |C i |-valued variables. In order to be multiplexed with symbols of S H , the lower priority bitstream b must be mapped into a sequence of |C i |-valued variables. Let us first consider the mapping of the realization s H = s 1 . . . s t . . . s KH into the sequence n 1 . . . n t . . . n KH . The sequence of n t -valued variables can also be seen as a unique Λ-valued variable γ, where Λ = KH t=1 n t . The variable γ hence verifies 0 ≤ γ < Λ. The quantity Λ denotes the number of different multiplexed sequences of codewords that can be used as a coded representation of the sequence s H . One of these sequences can be used as a multiplexed coded description of the two sequences s H and s L 1 . The sequence of multiplexed codewords to be transmitted will then depend on the bitstream representation b of the source S L .

The maximum amount of information that can be stored in the Λ-valued variable γ is theoretically log 2 (Λ) bits. However, since this variable is used to store bits of the bitstream b, only K ′ B = ⌊log 2 (Λ)⌋ bits can be stored, leading to an overhead equal to log2(Λ)

K ′ B . The last 2 K ′
B bits of b are then seen as the binary representation of the integer γ comprised between 0 and 2 K ′ B -1, that can be expressed as

γ = K ′ B r=1 b (r+KB -K ′ B ) 2 r-1 (2)
Note again that the variable γ is an index of the set of sequences that represent s H . The variable γ must then be expanded into a sequence of pairs (n t , q t ) that will provide entries in the multiplexed codes table. Let us recall that q t denotes the index of a codeword in the equivalence class associated to the realization s t of a given symbol S t . There are different ways to expand the variable γ into the sequence of n t -valued variables (q t ), where t = 1, . . . K H . One possible method is to use the recursive Euclidean decomposition: Definition 2: the recursive Euclidean decomposition of a positive integer x by a sequence y = (y 1 , y 2 , . . . , y K ) of 1 In fact, only a part of s L , whose length depends on the multiplexing capacity of the sequence s H 2 Any other subset of K ′ B bits of b can be used Algorithm 1 Conversion of the integer γ into the sequence q for t = 1 to K H do positive integers is the unique sequence r = (r 1 , r 2 , . . . , r K ) ∈ N K such that:

q t = γ ′ modulo n t γ ′ = γ ′ -
x = r 1 + y 1 (. . . (r t + y t (. . . (r K-1 + y K-1 r K ) . . .)) . . .)) ∀t ∈ [1..K], 0 ≤ r t ≤ y t -1 (3) It leads to expand γ as q 1 + n 1 (q 2 + n 2 (. . . (q KH -1 + n KH -1 q KH ) . . .)). Since each codeword index q t verifies q t = γ t-1 t ′ =1 n t ′ mod n t , (4)
the components q t can be calculated by Algorithm 1. In section III-C, we will see that these quantities can be computed with a lower complexity. In summary, the encoding proceeds as follows:

1) For t = 1, . . . K H , let n t be the cardinality of the equivalence class associated to the realization s t of the symbol S t .

2) The integer Λ is computed as Λ = KH t=1 n t . The number K ′ B of bits of the lower priority bitstream b that can be jointly encoded with the sequence S H is given by ⌊log 2 (Λ)⌋.

3) The K ′ B last bits of the low priority bitstream b are converted into the sequence of pairs (n t , q t ), t = 1 . . . K H , using Algorithm 1 of the Euclidean decomposition of (3). 4) The sequence of pairs (s t , q t ) provides the entries in the table of multiplexed codes, hence allows to select the sequence of c-bits fixed length codewords to be transmitted on the channel.

5) If K ′ B ≤ K B , the K B -K ′
B first bits of the bitstream b are then concatenated to the sequence of multiplexed codewords. Otherwise, the remaining symbols of S H are not multiplexed but sent using other codes on the channel, such as FLCs or Huffman codes. As all steps are reversible, the decoding proceeds in the reverse order. A schematic diagram of the encoding process, for the particular example proposed hereafter, is depicted in Fig. 1. 1. Link between the source S H , the low priority source formatted as a state flow q of |C i |-valued variables, and multiplexed codewords. Here, the bolded path denotes the chosen codewords.

a 3 a 3 a 1 a 2 s H a 1 a 5 a 2 a 4 (q t ) =26 b C 1 C 2 C 1 C 3 C 3 C 2 C 5 C 4 n 8 = 2 (n t ) s L Fig.
Example 2: let us consider the source alphabet described in section III and the multiplexed code given in Table I. Considering the sequence of symbols S H given in Table II, the algorithm proceeds first with the derivation of the sequence of n t variables as shown in Table II. This leads to Λ = log 2 (36) ≈ 2 5.17 . Therefore, the last 5 bits of the bitstream b = 11010 are used to process the variable γ = 5 r=1 b r 2 r-1 = 26. The Euclidean decomposition of the variable γ according to (1) leads to the sequence of q t variables (indexes of the codewords in the classes of equivalence associated to the sequence of symbols) given in Table II. The sequence of pairs (s t , q t ) provides directly the entries in the table of c-bits fixed length codewords, i.e. in the table of multiplexed codes.

C. Hierarchical decomposition of γ

The complexity order of the Euclidean decomposition of the global variable γ is O(K 2 H ) if Algorithm 1 is used. This overall complexity can be reduced by considering the following hierarchical decomposition of the variable γ. Note that this reduction in terms of complexity does not induce any overhead or redundancy. The resulting states q t are identical to the ones obtained by Algorithm 1. For sake of clarity, let us assume that K H = 2 l , with l ∈ N. The approach can be easily extended to any length. Let n t2 t1 denote the product t2 t=t1 n t . Similarly, q t2 t1 is defined as

q t2 t1 = ∆ q t1 + n t1 (q t1+1 + . . . (q t2-1 + n t2-1 q t2 ) . . .). (5)
If the number of terms in ( 5) is a power of two, i.e if t 2t 1 + 1 = 2 l ′ , then the entities n t2 t1 and q t2 t1 are respectively decomposed as

n t2 t1 = n t 2 +t 1 -1 2 t1 n t2 t 2 +t 1 +1 2 , (6)
q t2 t1 = q t 2 +t 1 -1 2 t1 + n t 2 +t 1 -1 2 t1 q t2 t 2 +t 1 +1 2 . ( 7
)
Algorithm 2 Conversion of the integer γ into the sequence of states q using the hierarchical algorithm (for sequences such that ∃ l/K H = 2 l ).

for j = 1 to l do {Processing of values n t2 t1 involved in further processing}

for i = 1 to 2 l-j do t 1 = (i -1) 2 j + 1 t 2 = i 2 j t m = (t 1 + t 2 -1)/2 n t2 t1 = n tm t1 n t2 tm+1 end for end for for j = l to 1 by -1 do {Processing of values q t2 t1 } for i = 1 to 2 l-j do t 1 = (i -1) 2 j + 1 t 2 = i 2 j t m = (t 1 + t 2 -1)/2 q tm t1 = q t2 t1 modulo n tm t1 q t2 tm+1 = (q t2 t1 -q tm t1 )/n tm t1

end for end for

Notice that γ = q 2 l 1 . Hence the state γ is decomposed as

γ = q 2 l-1 1 + n 2 l-1 1 q 2 l 2 l-1 +1 . (8)
Similarly, each term of ( 8) can be recursively decomposed using ( 6) and [START_REF] Maxted | Error recovery for variables length codes[END_REF]. Thus, the conversion of the integer γ into a sequence of states q can be done iteratively, as described in Algorithm 2. Note that each level j generates the variables required for level j + 1. The reverse algorithm follows naturally. It is shown in Appendix II that the complexity order of this algorithm is in O((K H ) r ), where r > 1 is the complexity exponent of the algorithms used for operations of large integers. Usually, the algorithm used is such that r goes from r = log3 log2 ≈ 1.58 to r = log9 log5 ≈ 1.37, depending on the length of the message [START_REF] Karatsuba | Multiplication of multidigit numbers on automata[END_REF] [START_REF] Zuras | More on squaring and multiplying large integers[END_REF]. Example 3: let us consider the decomposition proposed in Example 2. This decomposition can also be processed with the hierarchical algorithm as n 2 1 = 3, n 4 3 = 2, n 6 5 = 1, n 8 7 = 6, n 4 1 = 6, n 8 5 = 6 and q 4 1 = 26 mod n 4 1 = 2, q 8 5 = ⌊ 26

n 4 1 ⌋ = 4, q 2
1 = 2, q 4 3 = 0, q 6 5 = 0, q 8 7 = 8, q 1 = 1, q 2 = 0, q 4 3 = 0, q 3 = 0, q 4 = 0, q 5 = 0, q 6 = 0, q 7 = 1 and finally q 8 = 1.

D. Compression efficiency

Each symbol of the sequence S H is associated to a cbits fixed length codeword, leading to a total length of the multiplexed flow equal to c × K H . Therefore, the compression rate is determined by the amount of S L information that can be jointly coded by the sequence of multiplexed codewords. Let us consider again the realization s t of a symbol S t of the sequence S H . One can associate an n t -valued variable to the symbol value s t . The amount of data that can be represented by this variable is log 2 (n t ) bits. Since the c bits of the multiplexed codeword code both the symbol value s t and the index q t describing the low priority bitstream realization, the description length for the symbol s t is theoretically clog 2 (n t ) = -log 2 ( nt |X | ) bits. The amount of data, in bits, used to code the sequence s H is given by

- 1≤t≤KH log 2 ( n t |X | ). (9)
In order to minimize the mean description length (mdl) ĥ of the source S H , one has then to choose the partitioning C of the set of |X | fixed length codewords into equivalence classes

C i of cardinality |C i |, i = 1 . . . |A|, such that (|C 1 |, . . . |C i |, . . . |C |A| |) = argmin   - |A| i=1 µ i log 2 ( |C i | |X | )   = argmin( ĥ). ( 10
)
In order to have a rate close to the entropy bound of the source S H , the quantity |Ci| |X | should be as close as possible to µ i . We come back in section VI on how to calculate an efficient set of |C i | values for a given source pdf. Example 4: the entropy of the source S H introduced in section II is h = 2.122. The partition C given in Table I leads to an mdl of S H of ĥ = 2.166. Note that a Huffman encoder of the source would have led to an mdl equal to 2.2. The choice of c as well as of the partition C has an impact on the coding rate. For example, for c = 5, the partition into classes of respective cardinalities [START_REF] Takishima | Reversible variable length codes[END_REF][START_REF] Maxted | Error recovery for variables length codes[END_REF][START_REF] Capocelli | On the characterization of statistically synchronizable variable length codes[END_REF][START_REF] Bauer | Turbo fec/vlc decoding and its application to text compression[END_REF][START_REF] Bauer | Turbo fec/vlc decoding and its application to text compression[END_REF] produces an mdl equal to 2.124.

(|C 1 |, |C 2 |, |C 3 |, |C 4 |, |C 5 |) =
Property 1: ∀ ǫ > 0, ∃ c ∈ N such that ĥ -h ≤ ǫ
In other words, the mdl of S H can be made as close as required to the entropy, by increasing the length c of the codewords. The proof is provided in Appendix I.

IV. MULTIPLEXED CODES BASED ON A CONSTRAINED PARTITION C

The computational cost of the mapping (or multiplexing) algorithm can be further reduced by considering a constrained partition of the set of fixed length codewords. The preformatted lower priority bitstream is in this case not seen as the representation of a unique variable γ, but as a sequence of "elementary variables". An elementary variable is defined as a variable which takes its value in a set of small dimension: it is a k-valued variable, k being small.

A. Constrained choice of partition C

The |C i |-valued variables, 1 ≤ i ≤ |A|, are decomposed into a set of elementary variables as follows:

|C i | = νi j=1 f αi,j j , (11)
where f j are prime factors. The term f νi denotes the highest prime factor in the decomposition of |C i |. The term α i,j stands for the power of the prime factor f j in the decomposition

bit 1 bit 2 bit 3 3-valued variable 1 3-valued variable 2 0 0 0 0 0 0 0 1 0 1 0 1 0 0 2 0 1 1 1 0 1 0 0 1 1 1 0 1 1 2 1 1 0 2 0 1 1 1 2 1 2 2 TABLE III EXAMPLE OF TRANSFORMATION T (u T = 3).
of |C i |. This expansion of |C i |-valued variables can also be represented as:

|C i |-valued variable ⇄            α i,1 binary variables α i,2 3-valued variables . . . α i,νi f νi -valued variables (12)
In order to limit the encoding complexity, the partitioning of the set of c-bits fixed length codewords must be such that the set (|C i |) 1≤i≤|A| will not result in prime decompositions into factors greater than a given value f ν . These additional constraints on the partition C will induce an approximation of the pdf of the S H source. This approximation leads to a slightly higher mdl. We come back on this point in Section VIII-B.

B. Conversion of the low priority bitstream into f j -valued variables

Consecutive segments of the low priority bitstream are seen as the binary representations of integers to be decomposed this time into a set of binary, 3-valued, . . ., j-valued, . . ., f ν -valued variables. Let T be a transformation that takes u T bits and produces a set of f j -valued variables. Let v T ,j , 1 ≤ j ≤ ν, be the number of f j -valued variables produced by a transformation T . The number of possible values that can be taken by the set of f j -valued variables resulting from the transformation T applied on u T bits of the bitstream b is given by

Λ T = ν j=1 f vT ,j j . ( 13
)
This Λ T -valued variable allows to theoretically store log 2 (Λ T ) bits. Note that Λ T may be higher than 2 uT . Therefore, there is an overhead per bit given by

o T = log 2 (Λ T ) u T -1. ( 14
)
For example, the transformation T 15 shown in table III applied on 3 bits produces two 3-valued variables, leading to a number of states increased from 8 to 9.

Table IV enumerates the transformations used for f ν = 5. If f ν = 3, only the transformations that do not use 5-valued variables are valid: T 0 , T 3 , T 12 , T 15 , T 17 . They are sorted by increasing loss in compression efficiency. Among all the possible sets v T = (v T ,1 , . . . v T ,j , . . . v T ,ν ) explored under Algorithm 3 Choice of the transformations to be used z = 0 while sum(d j ) > 0 do {while some f j -valued variables are not computed yet} g Tz = ⌊min( dj vT z ,j )⌋ {Calculation of how many transformations T z can be used} for ∀j between 1 and ν do {Update the number of f jvariables to be transformed} d j = d j -g Tz * v Tz,j z = z + 1 {Try next transformation} end for end while previous constraints, only those for which the transformations introduce an overhead lower than 1% are kept. Increasing f ν reduces this overhead.

C. Optimal set of transformations

Let again s h be the realization of a source S H mapped into a sequence of n t -valued variables. Each n t -valued variable is decomposed into a set of f j -valued variables, each one being used α t,j times. Let d j be the total number of f j -valued variables involved in the expansion of the entire sequence of n t -valued variables. Let g Tz be the number of transformations T z , z = 0, . . . z max necessary to convert the low priority bitstream into the f j -valued variables. There are several possible sets g = (g T1 , . . . g Tz , . . . g Tz max ). The optimum set of transformations is the one which will minimize the total overhead, given by

O = zmax z=0 g Tz u Tz o Tz . ( 15
)
The choice of g is an optimization problem. Knowing the transformation set (T z ) 0≤z≤zmax and d = (d 1 , . . . d j , . . . d ν ), the optimal vector g is estimated by Algorithm 3.

t 1 2 3 4 5 6 7 8 st a 1 a 4 a 5 a 2 a 3 a 3 a 1 a 2 nt 3 1 1 2 1 1 3 2 α t,1 0 0 0 1 0 0 0 1 d 1 = 2 α t,2 1 0 0 0 0 0 1 0 d 2 = 2

D. Encoding procedure

Assuming that the low priority source s L has been preencoded using efficient VLCs, the encoding of the two sequences s H and s L using fast multiplexed codes proceeds as follows:

1) For each prime integer f j , 1 ≤ j ≤ ν, the number d j of f j -valued variables resulting from the expansion of the sequence of n t -valued variables associated to s H is calculated.

2) The number g Tz of transformations T z needed to convert the low priority bitstream b into a sequence of f j -valued variables is searched with the procedure described in Section IV-C. The transformations T z are then applied on consecutive segments of u Tz bits of the bitstream b, leading to the generation of sequences of f j -valued variables.

3) The n t -valued variable associated to the symbol realization s t has to be in bijection with several elementary f j -valued variables generated in the previous step. Then, at this point, the encoder has to choose the value q t taken by this n t -valued variable, such that the Euclidean decomposition of q t leads to the given sequence of f j -valued variables. This is done using the Euclidean multiplication. 4) The resulting pair (s t , q t ) provides entries in the multiplexed codewords table to select the codeword to be transmitted.

The decoding algorithm is made in the reverse order. As a consequence, it is not necessary to receive the whole sequence of codewords to decode the information of high priority s H : it is directly read in the multiplexed codes tables. Moreover, random access is also possible for the flow s H .

Example 5: let us consider the sequence of 8 high priority symbols given in Table V, with the corresponding mapping into the sequence of n t -valued variables according to the multiplexed codes given in Table I. The sequence of n t -valued variables is decomposed into a sequence of d 1 binary and d 2 3-valued variables with respective powers α t,1 and α t,2 (cf. table V). In a second step, one has to convert the low priority bitstream into the binary and 3-valued variables. It appears that both transformations T 0 and T 15 have to be used (see section IV-C) 2 times and 1 time, respectively. Thus, 5 bits can be multiplexed with the realization s H of S H . The transformation T 0 being identity, the first 3 bits of the bitstream b are multiplexed unchanged. For each transformation T , the u T input bits are seen as the binary representation of an integer γ T as

T Input bits → γ T → binary 3-valued T 0 : 1 → 1 → 1 T 0 : 1 → 1 → 1 T 15 : 010 → 2 → 02
and leading to the sequences 11 of binary variables and 02 of 3-valued variables. The value q t taken by the n t -valued variable associated to the symbol realization s t and to the realization of the low priority bitstream b is obtained from the sequences of binary (j = 1) and 3-valued (j = 2) variables. In the example above, the state q t is generated using either a bit for equivalence class C 2 , either a 3-valued variable for equivalence class C 1 . The resulting pair (s t , q t ) provides entries in the multiplexed codes table (see Table I), leading to the selection of multiplexed codes given in Table VI.

V. BINARY MULTIPLEXED CODES

The partition of the FLC into several equivalence classes can be chosen such that the storage capacity is an integer number of bits. It means that, for each class C i , a bijection exists between the interval of indexes [0, |C i |-1] and a varying integer number of bits, hence the following definition:

Definition 3: a binary multiplexed code is a multiplexed code such that ∀i ∈ [1..|A|], log 2 (|C i |) ∈ N.
Then, to each symbol a i of the source S H , one can associate an equivalence class C i , hence a set of fixed length codewords, which is in turn indexed by a sequence The number of bits that can be multiplexed with each symbol realization of the sequence S H is processed, which leads to segment b into (u 1 , . . . u t , . . . u KH ) = (1, ∅, ∅, 1, 0, 1, 0, 0), where notation ∅ indicates that the corresponding equivalence class is of cardinality 1 and is consequently not indexed by any bit. Then, for each high priority symbol index t, the couple (s t , u t ) indexes a codeword in the multiplexed table. The resulting multiplexed bitstream is 010 111 000 100 110 101 001 011.

U i = U 1 i . . . U
Instead of designing a FLC codebook for the high priority source S H , one can alternatively consider a variable length codetree. This VLC codetree can then be completed to form a binary FLC codetree (see Fig. 2). The symbols of the alphabet of the source S H will then be associated either to leaves of the FLC or to intermediate nodes. When they correspond to a node, all the corresponding leaves will constitute the equivalence class associated to the given symbol of S H . The cardinalities of the equivalence classes are the integer powers of 2. Hence, one can create a multiplexed code such that codewords of the same equivalence class have the same prefix. The suffix U i of the codewords can then be used to "store" the bits of the bitstream b (see Fig. 2).

Definition 4: a Binary multiplexed code derived from a VLC is a multiplexed code constructed such that (1) every prefix of a codeword is the representation of a symbol in a VLC tree, (2) every suffix corresponds to the multiplexed data of the low priority bitstream.

Example 7: considering the same sequences as in example 6, the VLC-based multiplexed code depicted in Fig. 2 The compression efficiency for the source S H is obviously given by the efficiency of the VLC considered. Similarly, since the low priority source before multiplexing has been preencoded, the corresponding compression efficiency depends on the VLC code considered for this source. The low priority source (e.g high frequencies of a wavelet decomposition) can be pre-encoded with an arithmetic coder.

VI. PROBABILITY DISTRIBUTION FUNCTION

APPROXIMATION

As already mentioned above, the mdl ĥ of the source S H is function of the sizes |C i | of the equivalence classes C i , i = 1 . . . |A|. This section describes an algorithm that aims at approximating, under constraints, the pdf of the source S H . The probability distribution function estimator μ is defined by

μ = (μ 1 , . . . μi , . . . μ|A| ) = |C 1 | |X | , . . . |C i | |X | , . . . |C |A| | |X | , (16) where
a i ∈ A m iff µ i < 1 |X | and a i ∈ A M iff µ i ≥ 1 |X | . ( 17
)
Let μ be the probability distribution function of the set of symbols a i ∈ A M . Each is given by

μi = µ i ai∈AM µ i (18)
The algorithm aiming at the best partition of the set of codewords proceeds as follows:

1) For each

a i ∈ A m , |C i | is set to 1.
2) The probability density function μ of symbols belonging to A M is calculated. Since |A m | codewords have already been affected to equivalence classes in the first step, there is still |X | -|A m | codewords to be assigned to equivalence classes. The expression [START_REF] Capocelli | On the construction of statistically synchronizable variable length codes[END_REF] can then be written as

(|C 1 |, . . . |C i |, . . . |C |A| |) = argmin log 2 (|X |) ai∈Am µ i - ai∈AM µ i log 2 ( |C i | |X |
) .

(19) Since the first part of expression ( 19) is a constant and |C i | = 1 when a i ∈ A m , the optimization is performed for symbols of A M only, and is expressed as

(|C i |) i∈AM = argmin -ai∈AM µ i log 2 ( |Ci| |X | ) = argmin -ai∈AM µi i∈A M µi (log 2 ( |Ci| |X |-|Am| ) +log 2 ( |X |-|Am| |X | )) = argmin -ai∈AM μi log 2 ( |Ci| |X |-|Am| ) -ai∈AM μi log 2 ( |X |-|Am| |X | ) = argmin -ai∈AM μi (log 2 ( |Ci| |X |-|Am| ) . ( 20
)
As ai∈AM μi = 1, expression (20) can be seen as the expression of the mdl when the index i is constrained to be such that a i ∈ A M . Consequently, assuming that |Ci| |X |-|Am| can take any value of interval [0, 1], the optimum is given by:

∀a i ∈ A M , |C i | |X | -|A m | = μi . ( 21
)
In the following, the index i is supposed to be chosen such that

a i ∈ A M . 3) For each a i ∈ A M , |C i | is set to the highest value in η such that |C i | ≤ μi (|X | -|A m |)
. However, at this point, some cardinalities may be equal to zero. In that case, the corresponding symbols are assumed to belong to A m instead of A M and the algorithm goes back to step 1. The remaining number of codewords to be assigned to equivalence classes is then given by

σ = |X | - |A| i=1 |C i | ≥ 0. 4) As a consequence, it is possible to increase some values of |C i |.
For each index i between 1 and |A|, if the value η k+1 is used instead of η k , the mdl ĥ decreases of a positive value µ i log 2 (

η k+1 η k ).
The decreasing mdl is induced by the assignment of η k+1 -η k pairs (a i , |C i |) to classes of equivalence and, in average per codeword, is given by

γ i = µi log2( η k+1 η k
) η k+1 -η k . These pairs are then sorted by decreasing values of γ i . The way to proceed is to choose the pair (a i , |C i |) with the highest associated value γ i . For a given variable |C i | set in the previous steps, it is possible to re-set it to η k+1 instead of η k only if η k+1 -η k ≤ σ, i.e. if there is a sufficient number of codewords still to be assigned. If it is not the case, the pair (a i , |C i |) is ignored in the next iteration of the algorithm. The procedure continues with the treatment of the pair with the next value

γ i . If η k+1 -η k ≤ σ, then |C i | = η k+1 , the value σ is set to σ -η k+1 + η k ,
and the variable γ i is updated for this symbol.

VII. ERROR HANDLING ON A BINARY SYMMETRIC

CHANNEL

In this section, we analyze the impact of channel errors on the distortion of the source S H . The part of the flow b that is not multiplexed is supposed to be lost if an error occurs on the channel. We consider a binary symmetric channel (BSC) with a bit error rate p. The pdf µ of S H being known, it is possible to calculate the reconstruction accuracy in terms of Symbol Error Rate (SER) or of Mean Square Error (MSE). Since each state q t calculated from b can be seen as a uniform random variable, a symbol a i of A has the same probability to be encoded with any codeword of its equivalence class C i . Hence, for a given codeword x ∈ C i ,

P (x) = µ i |C i | . ( 22
)
The probability to receive the symbol y, if x has been emitted, is function of the Hamming distance H(x, y) between the codewords x and y, and is given by P (y/x) = p H(x,y) (1p) c-H(x,y) . The mean distortion induced by the channel noise is then given by

E(∆ SH ) = x∈X P (x) y∈X P (y/x) ∆(x, y), (23)
where ∆(x, y) denotes the distortion of the source S H induced by the reception of y instead of x. Let respectively a i and a i ′ be the symbols decoded from x and y and ∆(a i , a i ′ ) the corresponding distortion measure. E.g., ∆(a i , a i ′ ) may be, respectively, the Kronecker operator or the quadratic error, depending on the error-resilience measure (respectively the SER or the MSE). This leads to

E(∆ SH ) = ai∈A µ i a i ′ ∈A ∆(a i , a i ′ ) 1 |C i | x∈Ci y∈C i ′ P (y/x), P (a i ′ /ai) (24)
where P (a i ′ /a i ) is the probability to decode the symbol a i ′ if a i has been transmitted. Therefore, the choice of the partition C has a major impact on the final distortion. However, the number of partitions for a given set of |C i |-values is very high, and is given by

|X |! |A| i=1 |Ci|!
. For example, the number of partitions such that

|C 1 | = 3, |C 2 | = 2, |C 3 | = 1, |C 4 | = 1, |C 5 | = 1, is 3 360.
Example 8: let us consider again the example of the source alphabet given in section III and the corresponding multiplexed code of codeword length c = 6. Let us define the alphabet values as a 1 = 1, a 2 = 2, a 2 = 3, a 3 = 4, a 4 = 5. We consider the two following distortion measures: the SER (in that case ∆(a i , a i ′ ) = 1 if a i = a i ′ , 0 otherwise) and the quadratic error (∆(a i , a i ′ ) = (a i -a i ′ ) 2 ). For a BSC of BER equal to 0.01, the summation in [START_REF] Huffman | A method for the construction of minimum redundancy codes[END_REF] leads to a SER of 0.0279 and an MSE of 0.0814 for a multiplexed code when a lexicographic index assignment is used. If the index assignment is optimized using a simulated annealing algorithm (e.g, [START_REF] Farvardin | A study of vector quantization for noisy channels[END_REF]), then the distortion is 0.0243 for the SER. However, we did not observe any improvement for the MSE criteria using this method (for this source).

For multiplexed codes derived from a VLC, the analysis of the SER simplifies as follows. One can remark that a given symbol a i is in error if and only if the prefix that identifies the equivalence class is in error. Hence, the probability that the symbol a i is in error is given by p c-log2(|Ci|) . This leads to the following expression of the SER:

SER SH = 1 - ai∈A µ i (1 -p) c-log2(|Ci|) . (25)
For low BER, the approximation

(1 -p) n = 1 -n p + O(p 2 ) leads to SER SH = p ai∈A µ i (c -log 2 (|C i )|) + O(p 2 ), (26)
where we recognize the expression of the mdl of the VLC code from which the multiplexed code has been derived: SER ≈ p × mdl. This SER is lower than for FLCs (p × mdl ≤ p × c), but the error-resilience is not uniformly better: symbols with higher probabilities of occurrence have a lower probability to be erroneous whereas symbols with lower probabilities have a higher probability to be erroneous. Note on the error-resilience of the low priority bitstream: in the general case and from equation ( 3), it appears that any error occurring on a multiplexed codeword at symbol clock t engenders, at least, an error on the quantity n t or on the state q t , or on possibly both. Consequently, the remaining states q t+1 • • • q KH are corrupted. Thus, the conversion of S L into a sequence of state is a step that is sensitive to bit errors. As a first approximation, the expectation of the position of the first impaired bit on the BSC is given by (cmdl SH )/(c × p) (instead of 1/p if S L is not multiplexed).

In binary multiplexed codes, the bits of S L are directly indexed by codewords. It follows that either a bit error may lead to drop or to insert some bits, either the number of decoded bits is correct but these bits are erroneous. Consequently, the bitstream is not fully desynchronized in terms of Levenshtein distance. This fact has been validated by simulations results (see section IX).

VIII. DISCUSSION AND PARAMETERS CHOICE

A. Redundancy in terms of the ratio S H versus S H +S L

We have seen above that c must be high enough so that one can find a partition that will closely approximate the pdf of the source S H . On the other hand, a high value for the parameter c will result in a decreased rate of synchronous data. Thus a compromise has to be found between reaching the entropy (compression efficiency) and the rate of synchronous data in the multiplexed flow (error-resilience). Fig. 3 shows the overhead (or redundancy) in terms of the ratio between the source S H mdl and the total amount of transmitted bits (S H +S L ). It can be seen that, if we consider multiplexed codes for joint encoding of motion and texture information in a video stream, for which realistic ratios are around 15%, the redundancy of the code is negligible (around 10 -6 ).

Fig. 3 shows the performance of multiplexed codes versus Huffman codes when both codes are designed for product alphabets. One can observe in Fig. 3 that when considering a ratio of S H versus the total number of symbols of 60%, the overhead remains rather low (10 -3 ). To achieve the same overhead with a Huffman code, one has to take the symbols by groups of 6 symbols (alphabet A 6 ). Fig. 3 also shows that using product alphabets (here A 2 , A 4 and A 6 ) allows to find a better trade-off between compression efficiency and the ratio of the high priority source versus the overall number of symbols, however at the cost of larger codetrees or tables. Note that using product alphabets preserves the strict-sense synchronization of the multiplexed codes. However, this hard synchronization is in this case guaranteed every n symbols, where n is the power of the alphabet.

B. Pdf approximation inherent to codes with constrained partitions

Let us consider the codes based on a constrained partition presented in Section IV-A. The ability to closely approximate the pdf of the source S H (hence the closeness to entropy) depends also on the parameter ν. The choice of the parameters c and ν has hence a strong impact on the performance, in terms of compression and error-resilience, of the codes. Let h be the entropy per symbol of the source and ĥ the mean codeword length produced by the algorithm. The overhead rate is given by ĥ-h h . As c bits are always used to represent a symbol of the source S H , the rate of high priority data in the final multiplexed c-bits codewords stream is given by h c . Fig. 4 shows the influence of the parameters c and ν on the rate ĥ-h h . For a source of cardinality |A| = 256, the value f ν = 5 leads to a small overhead when the codeword length is c = 14. Notice that the number of available transformations and the number of flows of f j -valued variables increases with ν. Therefore, the computational cost of the algorithm increases with ν.

The problem of finding good cardinalities for binary multiplexed codes is equivalent to the problem of finding good VLC binary codetrees. Thus, choosing the codewords length is equivalent to minimizing [START_REF] Capocelli | On the construction of statistically synchronizable variable length codes[END_REF] [24]. However, in order to meet the requirements in terms of the ratio of synchronous data, one may limit the length of the longest codeword [START_REF] Voorhis | Constructing codes with bounded codeword lengths[END_REF], which amounts to limit the parameter c. For this purpose, the algorithm proposed in [START_REF] Larmore | A fast algorithm for optimal length-limited huffman codes[END_REF] may be used.

C. Elements of encoding and decoding complexity

The conversion of the low priority source into a sequence of states is the step that has the highest computing cost. For binary multiplexed codes, we have seen in section V that this conversion is straightforward. In that case, the subsequent computing costs required for encoding and decoding are very similar to the ones involved in encoding and decoding VLCs. The encoding and decoding of multiplexed codes based on constrained partitions (see Section IV) has a higher computing cost but the order of complexity remains the same: the number of operations involved in the conversion of the low priority source is linear as the length of the sequence to encode increases. In section III-C, we have shown that the complexity of the general case is sub-quadratic, but above linear complexity. The exact complexity depends on the algorithms used for large integer multiplications and divisions.

Fig. 5 depicts quantitative results of encoding/decoding complexity using an implementation of Algorithm 2 based on the GNU Multiple Precision Arithmetic Library [START_REF]The gmp library[END_REF]. The computing cost, expressed in seconds, is compared with a typical implementation of Huffman codes. It is shown that the complexity of multiplexed codes, as the sequence length increases, remains tractable for most applications. Multiplexed Codes complexity has also been compared against the complexity obtained with soft decoding of VLCs based on a bit/symbol trellis [START_REF] Bauer | Turbo fec/vlc decoding and its application to text compression[END_REF] and on the BCJR algorithm [START_REF] Bahl | Optimal decoding of linear codes for minimizing symbol error rate[END_REF]. Fig. 5 shows that this solution has a rather higher complexity than multiplexed codes. This results from the fact that the corre-sponding trellis has a number of states which is quadratic as the length of the sequence increases. Note that many authors have proposed methods to decrease the complexity of VLCs, using pruning techniques or suboptimal trellis, e.g. [START_REF] Weidmann | Reduced-complexity soft-in-soft-out decoding of variable length codes[END_REF], but then the estimation is not optimal.

IX. SIMULATION RESULTS

The multiplexed code are all the most beneficial when applied to sources of different levels of priority. They intrinsically support unequal error protection of the corresponding sources. This is going to be illustrated below with a simple image coding system.

However, these codes can also be applied in the case, there is a unique source. They allow to guarantee hard synchronization for at least part of the sequence of symbols to be transmitted. To illustrate the advantage in such a situation, the performance of the multiplexed codes referred to as MUX codes in the sequel and on the figures) is first assessed in terms of SER (Hamming distance) and Levenshtein distance by considering the simple theoretical source given in Example 1 [START_REF] Maxted | Error recovery for variables length codes[END_REF]. The performance of the MUX code is compared to the performance achieved with the code V LC = {01, 00, 11, 100, 101} referred to as C5 in [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF]. We have chosen to use this code as a basis for comparisons since, among the 16 optimal codes described in [START_REF] Zhou | Synchronization recovery of variable length codes[END_REF], this is the code which leads to the best performance in terms of Levenshtein distance. For the experiments reported for the theoretical source, we have considered sequences of 1000 symbols. The sequence of 1000 symbols is divided into two segments. The first and second segments are assumed to be the sources of high and low priority respectively. The second segment of the sequence of symbols is thus coded with this code referred to as V LC on the figures. The sequence is divided in such a way that the multiplexing capacity of the high priority source is greater or equal to the length of the bitstream representation of the low priority source. The multiplexed code used in this first experiment is the lexicographic binary multiplexed code such that |C

1 | = 4, |C 2 | = 1, |C 3 | = 1, |C 4 | = 1 and |C 5 | = 1.
This code is derived from the VLC code {0, 100, 101, 110, 111}.

This code leads to a ratio between the high priority symbols to the overall number of symbols equal to 0.72. The mdl measured is close to the theoretical value of 2.2 with both the MUX and the V LC (C5) codes. The results have been averaged over 20000 simulations. For both codes, we have respectively measured the overall performances in terms of normalized Levenshtein distance and SER. Fig. 6 and Fig. 7 show the performance in terms of normalized Levenshtein distance and SER obtained on the entire source with both codes (curves VLC:S; MUX:S). These figures also show the normalized Levenshtein distance and SER obtained on the high (MUX:S H ) and low (MUX:S L ) priority symbols separately. It can be observed that the Levenshtein distance values averaged on the entire sequence are very close with both codes. However, Fig. 7 shows the advantage of multiplexed codes if we consider strict-sense synchronization, i.e. the SER or Hamming distance. It can also be observed that significantly lower normalized Levenshtein distance and SER values are obtained for the symbols to be considered as high priority symbols. One can also verify on the curves that the SER approximation given in section VII (SER SH ≈ 2.2 × BER) is valid. The benefit of this inherent unequal error protection feature is better illustrated below with a concrete image coding example.

The multiplexed codes have then been experimented with real sources to show the benefits of both the inherent unequal error protection feature and of the hard synchronization property guaranteed for the high priority symbols. We have considered a simple image coding system where the image is first decomposed into 16 subbands using a two-stage wavelet transform. The low and high frequency subbands have been quantized respectively on 7 and 4 bits. The high frequencies have been encoded with an Huffman algorithm, and then multiplexed into the low frequencies using lexicographic multiplexed codes of parameters c = 16 and f ν = 5. A synchronization marker has been used every four lines for both Huffman codes and the low priority part of the multiplexed codes. The bit rates obtained with Huffman and multiplexed codes are quite similar and given respectively by 1.713 bit per pixel (bpp) and 1.712 bpp. When considering FLCs, the average bit rate obtained is 6.187 bpp. These FLCs are based on the lexicographic index assignment. Fig. 8 depicts the PSNR and visual qualities obtained with the three codes. This figure evidences significant improvements both in PSNR and visual quality when using multiplexed codes even in comparison with FLCs and for a similar compression factor as the one obtained with Huffman codes.

X. CONCLUSION

This paper describes a new family of codes called "multiplexed codes". These codes avoid the "de-synchronization" phenomenon for symbols considered to be of high priority, and allow thus to confine error propagation to symbols assumed to be of lower priority. The low priority symbols can be preencoded with any efficient VLC. A FLC is created for the high priority source, its inherent redundancy being exploited to represent information from the low priority stream. These codes are shown to reach asymptotically the entropy bound for both (low and high priority) sources. Several methods of multiplexing of the two sources are described. Sub-classes of multiplexed codes relying on constrained partitions of the FLC codebook are also introduced. The codes called binary multiplexed codes can be constructed from classical VLC. They thus inherit the compression properties of the VLC considered and allow for multiplexing with very low complexity. Theoretical and simulation results, using both theoretical and real sources, show that these codes are errorresilient at almost no cost in terms of compression efficiency. Another key advantage is that they allow to make use of simple decoding techniques. They hence appear to be excellent alternatives to reversible variable length codes (which suffer from some penalty in terms of compression efficiency, while not avoiding completely error propagation despite a decoding in two passes) or to classical VLCs for which robust decoding makes often use of computationally expensive Bayesian estimation techniques.

log2|Ci| i of log 2 (Fig. 2 .

 22 Fig. 2. Codetree describing both S H and B

|A|

  i=1 |C i | ≤ |X |. Let η = {η 1 = 1, . . . η k , . . .}be the set of possible values |C i |. A higher number of valid values for |C i | leads indeed to a better approximation of the pdf of the source S H . The probability µ i of a symbol cannot be lower than 1 |X | . The alphabet A is partitioned into two subsets A m and A M , defined respectively by

Fig. 4 .

 4 Fig. 4. Closeness to entropy of fast-computable Multiplexed codes (Gaussian source).

Fig. 6 .Fig. 7 .

 67 Fig.6.Performances in terms of normalized Levenshtein distance of multiplexed codes vs Huffman codes.

  .50 dB, bpp = 6.187 PSNR = 24.78 dB, bpp = 1.713 PSNR = 30.71 dB, bpp = 1.712 PSNR = 16.79 dB PSNR = 15.59 dB PSNR = 24.82 dB PSNR = 13.99 dB PSNR = 12.64 dB PSNR = 19.63 dB

Fig. 8 .

 8 Fig. 8. PSNR performance and visual quality obtained respectively with FLCs, Huffman codes and multiplexed codes. The channel bit error rates are 0.0005 (top images), 0.005 (middle images) and 0.05 (bottom images).

qt nt end for

	t	st	nt	γ	qt	codeword
	1	a 1	3	26	2	001
	2	a 4	1	8	0	000
	3	a 5	1	8	0	111
	4	a 2	2	8	0	011
	5	a 3	1	4	0	110
	6	a 3	1	4	0	111
	7	a 1	3	4	1	000
	8	a 2	2	1	1	001
				TABLE II
	EUCLIDEAN DECOMPOSITION OF THE VARIABLE γ AND CORRESPONDING
		MULTIPLEXED CODEWORDS.

TABLE V CALCULATION

 V OF THE NUMBER OF AVAILABLE f j -VALUED VARIABLES.

  leads to the multiplexed bitstream below.

	s t	a 1	a 4	a 5	a 2	a 3	a 3	a 1	a 2
	prefix	00. 110 111 01. 10. 10. 00. 01.
	u t	..1	...	...	..1 ..0 ..1 ..0 ..0
	codeword 001 110 111 010 100 101 000 010

  . The source is the one given in Example 1. Multiplexed codes designed for the product alphabets {a 1 , a 2 , a 3 , a 4 , a 5 } n with n = 1, 2, 4, 6 have been considered. Compression efficiencies of corresponding Generalized Huffman codes are also depicted.

		1		MUX: A 6 MUX: A 4 MUX: A 2 MUX: A
		0.1
		0.01		Hu man: A Hu man: A 2 Hu man: A 4			c=3
	overhead	0.001		Hu man: A 6	c=6
		0.0001
		1e-05
		1e-06	0	0.2	0.4 rate of SH	0.6	0.8	1
	Fig. 3. Compression efficiency vs Ratio	S H S H +S L

  by assigning the cardinalities|C i | , under the constraints log 2 (|C i |) ∈ N and ai∈A |C i | = 2 c. It is well known that the algorithm that will provide the best approximation of the source S H pdf is the Huffman algorithm

		10
		1
		0.1
	Processing time (in seconds)	0.001 0.01
		0.0001
							MUX encoding
							MUX decoding
		1e-05				Huffman encoding
						Huffman decoding
						Huffman soft decoding
		1e-06	10	100	1000	10000	100000	1e+06
					Sequence length
	Fig. 5. Encoding and decoding processing durations of Multiplexed codes
	(c=3) based on Algorithm 2. Source of Example 1 has still been considered.
	Corresponding durations are provided for encoding and decoding (hard and
	soft) of Huffman codes.

XI. ACKNOWLEDGMENT

The authors would like to thank the anonymous reviewers for their very constructive and helpful comments.

APPENDIX I PROOF OF PROPERTY 1

Let ǫ > 0. The problem consists in finding a partition C that verifies the property. Let c h be the integer defined as c h = ⌈-log 2 (min i∈A (µ i ))⌉. By definition,

For all c ∈ N such that c ≥ c h , we choose

The difference δ h between the mdl and the entropy is given by

Therefore, from ( 28), it can be seen easily that δ h verifies

From equations ( 27) and (30) we get

Notice that, for any code for which the inequality c ≥ c h is true, this proof provides an upper bound for the mdl.

APPENDIX II COMPLEXITY OF THE HIERARCHICAL DECOMPOSITION OF THE VARIABLE γ

The complexity of the hierarchical decomposition of the variable γ strongly depends on the complexity of the algorithms used for the multiplication, division and modulo operations. Karatsuba [START_REF] Karatsuba | Multiplication of multidigit numbers on automata[END_REF] has shown that these operations have a subquadractic complexity. For example, the multiplication can be done in O(N log3 log2 ) with the Karatsuba algorithm. For implementation purpose, efficient algorithms are described in the GNU Multiple precision computing library [START_REF]The gmp library[END_REF]. Now, let us assume that the complexity of long integer operations are given by C(K), where K is the block size of the largest variable involved. Then, processing the level j in Algorithm 2 has the complexity 2 l-j C(2 j ) = K H 2 -j C(2 j ). Assuming that C(K) = O(K r ) with r > 1, the overall complexity can be written as

Hence, from C(K) = O(K r ), we deduce that ∃A > 0 ∃B > 0 / ∀K > 0, C(K) < A + B K r . This leads to the following inegality

< K H A + K H B 2 (l+1)(r-1) .

(

Since we have

it appears that the complexity of processing the level l and the overall complexity are of the same order, i.e O((K H ) r ).


Serena Nik-Zainal

Reiner Strick

Mekayla Storer

Ni Huang

Roland Rad

Lionel Willatt

Tomas Fitzgerald

Vicki Martin

Richard Sandford

Nigel P Carter

Andreas R Janecke

Stefan P Renner

Patricia G Oppelt

Peter Oppelt

Christine Schulze

Matthew Hurles

Matthias W Beckmann

Pamela L Strissel

Charles Shaw-Smith

High incidence of recurrent copy number variants in patients with isolated and syndromic Mullerian aplasia

Keywords: Copy number variation, Mullerian aplasia, MURCS, MRKH

Background: Congenital malformations involving the Mullerian ducts are observed in around 5% of infertile women. Complete aplasia of the uterus, cervix, and upper vagina, also termed Mullerian aplasia or Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome occurs with an incidence of around 1 in 4,500 female births, and occurs in both isolated and syndromic forms.

Previous reports have suggested that a proportion of cases, especially syndromic cases, are caused by variation in copy number at different genomic loci.

Methods: In order to obtain an overview of the contribution of copy number variation to both isolated and syndromic forms of Mullerian aplasia, we performed copy number assays in a series of 63 cases, of which 25 were syndromic and 38 isolated.

Results: We report a high incidence (9/63, 14%) of recurrent copy number variants in this cohort. These comprised four cases of microdeletion at 16p11.2, an autism susceptibility locus not previously associated with Mullerian aplasia, four cases of microdeletion at 17q12, and one case of a distal 22q11.2 microdeletion. Microdeletions at 16p11.2 and 17q12 were found in 4/38 (10.5%) cases with isolated Mullerian aplasia, and at 16p11.2, 17q12 and 22q11.2 (distal) in 5/25 cases (20%) with syndromic Mullerian aplasia.

Conclusion:

Our finding of microdeletion at 16p11.2 in 2/38 (5%) of isolated and 2/25 (8%) of syndromic cases suggests a significant contribution of this copy number variant alone to the pathogenesis of Mullerian aplasia. Overall, the high incidence of recurrent copy number variants in all forms of Mullerian aplasia has implications for our understanding of the aetiopathogenesis of the condition, and for genetic counseling in families affected by it.

Introduction

The incidence of congenital malformations involving the Mullerian ducts in the general population is ~5 per 1000, and of infertile women more frequent at 35-63 per 1000 [START_REF] Raga | Reproductive impact of congenital Mullerian anomalies[END_REF][START_REF] Nahum | Uterine anomalies. How common are they, and what is their distribution among subtypes?[END_REF].

Among the most common uterine anomalies are uterine duplications, uterine indentations, and partial uterine aplasias, which occur as a result of incomplete Müllerian and Wollfian duct fusions during development. The Mayer-Rokitanski-Küster-Hauser (MRKH) syndrome (OMIM: 277000), occurring in 1 of 4500 live female births, is the most severe type where complete aplasia of the uterus, cervix, and upper vagina is found leading to failure to menstruate and infertility, despite normal secondary sexual characteristics [START_REF] Buttram | Mullerian anomalies: a proposed classification. (An analysis of 144 cases)[END_REF][START_REF] Folch | Mullerian agenesis: etiology, diagnosis, and management[END_REF][START_REF] Morcel | Programme de Recherches sur les Aplasies Mulleriennes, Guerrier D. Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome[END_REF].

The different MRKH subtypes are clinically classified into type I (typical or isolated) patients with normally developed fallopian tubes, ovaries and urinary tract, and type II (atypical) patients, with Fallopian or ovarian abnormalities, and additional malformations, which typically involve the urinary tract and spine [START_REF] Strubbe | Mayer-Rokitansky-Kuster-Hauser syndrome: distinction between two forms based on excretory urographic, sonographic, and laparoscopic findings[END_REF][START_REF] Oppelt | Clinical aspects of Mayer-Rokitansky-Kuester-Hauser syndrome: recommendations for clinical diagnosis and staging[END_REF]. The acronym MURCS (Müllerian-Renal-Cervicothoracic Somite Abnormalities, OMIM: 601076) applies to some of these cases. Craniofacial (dysmorphism, microtia) and cardiovascular malformations, and learning difficulties/mental retardation may also be associated [START_REF] Oppelt | Female genital malformations and their associated abnormalities[END_REF]. The reproductive and psychosocial consequences of the disorder are severe; despite this, from the aetiological perspective, it has been relatively poorly studied. Nonetheless, over the last decade, evidence has accumulated to suggest that genetic factors may be important.

Notwithstanding that patients with Mullerian aplasia are usually precluded from reproducing, familial clustering of the disorder has been reported, with apparently autosomal dominant inheritance [START_REF] Tiker | Familial mullerian agenesis[END_REF][START_REF] Shokeir | Aplasia of the Mullerian system: evidence for probable sex-limited autosomal dominant inheritance[END_REF]. For example, figures 6 and 9 (reference [START_REF] Shokeir | Aplasia of the Mullerian system: evidence for probable sex-limited autosomal dominant inheritance[END_REF]) show Mullerian aplasia segregating in apparently autosomal dominant fashion in two or more generations. In a small proportion of cases, a specific genetic aetiology has been identified. Mullerian aplasia with hyperandrogenism and renal malformations (OMIM 158330) is due to mutations in WNT4 [START_REF] Biason-Lauber | A WNT4 mutation associated with Mullerianduct regression and virilization in a 46,XX woman[END_REF].

Beside single gene involvement, there is evidence implicating copy number variation in the pathogenesis of Mullerian aplasia. Case reports in the literature have linked MURCS to microdeletions at 17q12 [START_REF] Bernardini | Recurrent microdeletion at 17q12 as a cause of Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome: two case reports[END_REF] and 22q11.2 [START_REF] Sundaram | Primary amenorrhea and absent uterus in the 22q11.2 deletion syndrome[END_REF]; in one small series of 14 cases, both of these loci were identified and two more suggested, a microduplication at 1q21.1 and microdeletion at Xq21.31 [START_REF] Cheroki | Genomic imbalances associated with mullerian aplasia[END_REF]. Mullerian aplasia has been reported in association with Thrombocytopenia-Absent Radius syndrome, due to microdeletion at chromosome 1q21.1 [START_REF] Klopocki | Complex inheritance pattern resembling autosomal recessive inheritance involving a microdeletion in thrombocytopenia-absent radius syndrome[END_REF].

Typically, syndromes due to recurrent copy number variants exhibit wide phenotypic variability, a fact which has been recognized since the early descriptions of the 22q11.2 microdeletion syndrome [START_REF] Ryan | Spectrum of clinical features associated with interstitial chromosome 22q11 deletions: a European collaborative study[END_REF] and which continues to be recognized in the new syndromes which have been described since the advent of high-resolution array-based studies (reviewed in ref [START_REF] Mefford | Duplication hotspots, rare genomic disorders, and common disease[END_REF]. A further example of this variability is given by the 16p11.2 microdeletion syndrome. This was initially described in cohorts of patients with autism [START_REF] Kumar | Recurrent 16p11.2 microdeletions in autism[END_REF][START_REF] Weiss | Association between microdeletion and microduplication at 16p11.2 and autism[END_REF]. Later, the recognition was made that patients with 16p11.2 microdeletions have a more complex phenotype than those patients with presumed multifactorial forms of autism spectrum disorder, including dysmorphism, congenital anomalies, growth disturbance, motor delay, and epilepsy [START_REF] Fernandez | Phenotypic spectrum associated with de novo and inherited deletions and duplications at 16p11.2 in individuals ascertained for diagnosis of autism spectrum disorder[END_REF][START_REF] Shinawi | Recurrent reciprocal 16p11.2 rearrangements associated with global developmental delay, behavioural problems, dysmorphism, epilepsy, and abnormal head size[END_REF]. Most recently, a strong association between 16p11.2 microdeletions and obesity was reported, especially where cognitive disability was also present [START_REF] Walters | A new highly penetrant form of obesity due to deletions on chromosome 16p11.2[END_REF]. A detailed explanation of how such marked phenotypic variability can arise from copy number variation at a single locus is lacking; environmental factors, epigenetic changes, and 'second hits' [START_REF] Girirajan | A recurrent 16p12.1 microdeletion supports a two-hit model for severe developmental delay[END_REF] may all play a part. For the present, it is likely that the full range of phenotypic variability of these disorders has yet to be explored, and that more associations may be identified through the study of different patient cohorts. Here, we report copy number analysis of DNA samples from a cohort of 63 individuals with Mullerian aplasia, of which 38 were classified as isolated or typical Mullerian aplasia (60 3%) and 25 were classified as syndromic or atypical Mullerian aplasia (39 7%), from which 11 had a diagnosis of MURCS association.

Results

Our study demonstrated a strikingly high incidence of recurrent copy number variants, with 9 (14%) of the 63 samples studied having a copy number variant of this type (Table 1 and Figures 1 and2). These CNVs were confirmed by quantitative PCR (Figure 3). For isolated Mullerian aplasia, this applied to four of 38 (10 5%) of patients, and, for syndromic Mullerian aplasia, to five of 25 (20%). The results are summarized in Tables 123. In isolated Mullerian aplasia, we identified two microdeletions at 17q12, and two at 16p11.2. Of the five syndromic cases, two had microdeletions at 17q12, two had microdeletions at 16p11.2, and one had a 'distal' 22q11.2 microdeletion. A sixth case had a microduplication at 2q11.2, a locus recently reported to represent a novel recurrent copy number variant, though as yet phenotypic data for this new disorder have not been provided [START_REF] Rudd | Segmental duplications mediate novel, clinically relevant chromosome rearrangements[END_REF]. Interpretation in this case is made additionally difficult by the presence of a large (4 6 Mb) deletion at 2p24.3. Either, or conceivably both, of these imbalances may have contributed to the phenotype in this case. There were two instances of microdeletion at 16p11.2 in the control cohort; no instances of microdeletion at 17q12 and 22q11.2 (distal) were identified.

We observed that 4/63 (6 %) of cases in our series had microdeletions at 16p11.2, the first time that this locus has been associated with Mullerian aplasia. Recent reports have identified microdeletions at this locus in cohorts of patients with autism spectrum disorder [START_REF] Weiss | Association between microdeletion and microduplication at 16p11.2 and autism[END_REF] and obesity [START_REF] Walters | A new highly penetrant form of obesity due to deletions on chromosome 16p11.2[END_REF]. These patients often have syndromic features such as dysmorphic facial features and some congenital malformations including vertebral anomalies, [START_REF] Fernandez | Phenotypic spectrum associated with de novo and inherited deletions and duplications at 16p11.2 in individuals ascertained for diagnosis of autism spectrum disorder[END_REF][START_REF] Shimojima | A familial 593-kb microdeletion of 16p11.2 associated with mental retardation and hemivertebrae[END_REF] and in one case, micropenis [START_REF] Fernandez | Phenotypic spectrum associated with de novo and inherited deletions and duplications at 16p11.2 in individuals ascertained for diagnosis of autism spectrum disorder[END_REF], but not female reproductive tract malformations. The enrichment of this locus in patients with Mullerian aplasia compared with controls (4/63 patients compared with 2/7,366 controls) is highly statistically significant (p = 6 96e-8, Fisher's exact test) and suggests a hitherto unappreciated role for genes within the deletion interval in the development of the Mullerian derivatives.

In keeping with previously published studies [START_REF] Bernardini | Recurrent microdeletion at 17q12 as a cause of Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome: two case reports[END_REF][START_REF] Cheroki | Genomic imbalances associated with mullerian aplasia[END_REF] we found microdeletions at 17q12 in patients with both syndromic and apparently isolated Mullerian aplasia. This finding is in line with previous reports in the literature linking this locus to syndromic Mullerian aplasia [START_REF] Cheroki | Genomic imbalances associated with mullerian aplasia[END_REF] and, in a single case, to apparently isolated Mullerian aplasia [START_REF] Bernardini | Recurrent microdeletion at 17q12 as a cause of Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome: two case reports[END_REF].

Several case reports have demonstrated an association between the velocardiofacial syndrome (VCFS) associated microdeletions at 22q11.2 and Mullerian aplasia [START_REF] Cheroki | Genomic imbalances associated with mullerian aplasia[END_REF][START_REF] Sundaram | Primary amenorrhea and absent uterus in the 22q11.2 deletion syndrome[END_REF][START_REF] Devriendt | Chromosome 22q11 deletion presenting as the Potter sequence[END_REF]. Recently, a novel genomic disorder was reported due to microdeletions at an adjacent, telomeric, locus, and this disorder was given the name '22q11.2 distal deletion' [START_REF] Ben-Shachar | 22q11.2 distal deletion: a recurrent genomic disorder distinct from DiGeorge syndrome and velocardiofacial syndrome[END_REF]. The complex genomic architecture at this locus gives rise to 'nested' microdeletions within the critical interval. In the original report of this syndrome [START_REF] Ben-Shachar | 22q11.2 distal deletion: a recurrent genomic disorder distinct from DiGeorge syndrome and velocardiofacial syndrome[END_REF], attention was drawn to two cases with distal nested microdeletions, one of which had an isolated congenital heart defect. We now describe a case with syndromic Mullerian aplasia/MURCS and the same, distal nested microdeletion.

We identified a microduplication at 2q11.2 in a patient with features of MURCS association.

Copy number variation at this locus has recently been reported and the suggestion made that this constitutes a novel genomic disorder [START_REF] Rudd | Segmental duplications mediate novel, clinically relevant chromosome rearrangements[END_REF]. However no phenotypic data concerning this duplication are available, and the interpretation in our case is further confounded by the co- The intervals delineated by these copy number variants harbour in some cases interesting candidate genes. The approximately 0 55 Mb interval at 16p11.2 contains TBX6, a gene previously implicated in development of paraxial mesoderm [START_REF] White | Defective somite patterning in mouse embryos with reduced levels of Tbx6[END_REF] but with no known role in formation of the Mullerian ducts. There are no other compelling developmental candidate genes in the critical interval. The 1.4 Mb interval at 17q12 harbours two genes with known roles in reproductive tract development: HNF1B, in which mutations have been described in patients with renal cysts and diabetes [START_REF] Lindner | A novel syndrome of diabetes mellitus, renal dysfunction and genital malformation associated with a partial deletion of the pseudo-POU domain of hepatocyte nuclear factor-1beta[END_REF]; these patients also have genital tract malformations such as bicornuate uterus and uterus didelphys, but absence of uterus and fallopian tubes has been reported [START_REF] Lindner | A novel syndrome of diabetes mellitus, renal dysfunction and genital malformation associated with a partial deletion of the pseudo-POU domain of hepatocyte nuclear factor-1beta[END_REF]. Secondly, this region harbours LHX1. Mutations in this gene have not been described in humans, but mice with targeted knockout of Lhx1 have absence of uterus and oviducts [START_REF] Kobayashi | Requirement of Lim1 for female reproductive tract development[END_REF]. Sequencing of LHX1 in patients with Mullerian aplasia has to date not revealed any mutations (ref [START_REF] Bernardini | Recurrent microdeletion at 17q12 as a cause of Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome: two case reports[END_REF], and our unpublished data). The 'distal 22q11.2' locus harbours just four genes, RTDR1, RAB36, GNAZ and BCR. Germline mutations in none of these four genes have been described in humans, and developmental malformations have not been reported in BCRnull [START_REF] Voncken | Increased neutrophil respiratory burst in bcr-null mutants[END_REF] or GNAZ-null mice [START_REF] Hendry | Hypertolerance to morphine in G(z alpha)-deficient mice[END_REF], for which mouse data are available.

11 non-genomic disorder type copy number variants were identified in the case cohort which were absent in controls, but none of these occurred in more than one patient, and the absence of parental samples is an additional factor limiting our interpretation of their significance. Table 3 lists copy number variants occurring in a single patient in the case cohort and not in the control cohort. Three of these, at 2p24.1-24.3, 15q21.1 and 18q23, were selected for validation by qPCR; all three were confirmed (Figure 3). The deletion of 4 6 Mb at 2p24.1-24.3 occurred in a patient with a double segment imbalance, the other copy number variant being a microduplication at 2q11.2, discussed above. One imbalance harbouring a potentially interesting candidate gene was noted, a 200 Kb duplication at 18q23 encompassing SALL3, a member of the SALL gene family, of as yet unknown function. Replication of this study on larger cohorts will be needed in order to determine whether these single occurrence CNVs are significant.

Discussion

Previous reports have identified copy number variants in Mullerian aplasia [START_REF] Bernardini | Recurrent microdeletion at 17q12 as a cause of Mayer-Rokitansky-Kuster-Hauser (MRKH) syndrome: two case reports[END_REF][START_REF] Cheroki | Genomic imbalances associated with mullerian aplasia[END_REF], but these have been case reports and small series, and studies of medium or large scale series have to date not been carried out. Here, we give results of copy number analysis of a series of 63 patients, identifying a strikingly high incidence of 9/63 (14%) of previously characterized microdeletions. Microdeletions at 16p11.2 and 17q12 in particular are highly enriched in the case population in comparison to the control population. Both of these microdeletions, and the 22q11.2 distal microdeletion, have previously been associated with congenital malformations: of the spine in the case of 16p11.2 [START_REF] Shimojima | A familial 593-kb microdeletion of 16p11.2 associated with mental retardation and hemivertebrae[END_REF], the genito-urinary tract in the case of 17q12 [START_REF] Mefford | Recurrent reciprocal genomic rearrangements of 17q12 are associated with renal disease, diabetes, and epilepsy[END_REF] and cardiovascular system in the case of the distal 22q11.2 deletion [START_REF] Ben-Shachar | 22q11.2 distal deletion: a recurrent genomic disorder distinct from DiGeorge syndrome and velocardiofacial syndrome[END_REF]. Thus, although data on inheritance for these cases is lacking, there can be little doubt that these microdeletions are contributing to the phenotypes of isolated and syndromic Mullerian aplasia.

The apparent contribution of genomic disorder-type copy number variants to congenital malformations occurring in isolation is of particular interest. A recent study of non-syndromic tetralogy of Fallot, a complex congenital heart malformation, gives weight to the idea that genomic disorders may be associated with isolated congenital malformations [START_REF] Greenway | De novo copy number variants identify new genes and loci in isolated sporadic tetralogy of Fallot[END_REF]. Copy number variants corresponding to known genomic disorders were identified at 22q11.2 (two, both loss), and at 1q21.1 (four gain, one loss) in 7/512 (1.4%) of individuals. In our study, copy number variants associated with genomic disorders were identified in 4/43 (9%) of isolated Mullerian aplasia cases, a significantly higher figure than for isolated tetralogy of Fallot.

Our study extends the phenotypic variability associated with microdeletions at 16p11.2. The previously known phenotypic consequences of microdeletions at this locus are: autism spectrum disorder [START_REF] Kumar | Recurrent 16p11.2 microdeletions in autism[END_REF]; epilepsy [START_REF] Shinawi | Recurrent reciprocal 16p11.2 rearrangements associated with global developmental delay, behavioural problems, dysmorphism, epilepsy, and abnormal head size[END_REF]; developmental delay/learning disability and dysmorphism/congenital anomalies [START_REF] Fernandez | Phenotypic spectrum associated with de novo and inherited deletions and duplications at 16p11.2 in individuals ascertained for diagnosis of autism spectrum disorder[END_REF][START_REF] Shinawi | Recurrent reciprocal 16p11.2 rearrangements associated with global developmental delay, behavioural problems, dysmorphism, epilepsy, and abnormal head size[END_REF]; and obesity [START_REF] Walters | A new highly penetrant form of obesity due to deletions on chromosome 16p11.2[END_REF]. Both isolated and syndromic forms of Mullerian aplasia can now be added to this list. Increasingly, these findings raise the question of which factor or factors are responsible for determining the phenotypic outcome in patients with 16p11.2 microdeletions in particular and genomic disorders in general. It can be hypothesized that this microdeletion, and others like it, provides an early and general developmental perturbation, the ultimate and precise consequence of which is independent of the perturbation itself, but depends on other factors. Future studies may begin to address this question, through sequencing, identification of 'second hits' [START_REF] Girirajan | A recurrent 16p12.1 microdeletion supports a two-hit model for severe developmental delay[END_REF] and other modalities.

Our results have implications for genetic counseling in Mullerian aplasia. At present, few women with the isolated form of this malformation are referred to clinical geneticists.

Mullerian aplasia results in infertility, and so it might be argued that the need for genetic counseling is less (although it might be important in cases where egg donation is being considered). More critically, though, the diagnosis of a microdeletion at 16p11.2 , 17q12 or 22q11.2 (distal) has potential implications for other family members. Copy number variants associated with genomic disorders may be inherited from a phenotypically normal parent and transmitted to other family members. We anticipate that families would wish to be made aware of the associations with autism (16p11.2) or young-onset diabetes and renal cysts (17q12) and cardiac malformations (22q11.2 (distal)) as well as with Mullerian aplasia, notwithstanding that our poor understanding of the penetrance and variable expressivity of these disorders makes genetic counseling difficult. In conclusion, our data support the contention that detailed copy number assays should be carried out in the assessment of both isolated and syndromic forms of Mullerian aplasia. We feel that the high incidence of recurrent copy number variants in these patients makes it reasonable to recommend that they should be referred to a clinical geneticist for assessment, and, where appropriate, for genetic testing and family counselling.

Materials and Methods

Patient samples

Patient samples and phenotypic data were collected by the Department of Obstetrics and Gynecology, Erlangen, Germany. The project was approved by the Ethics Review Board of Friedrich-Alexander-University, Erlangen-Nuremberg, Germany. Informed consent for genetic studies was obtained in each case.

Copy number assay

Copy number analysis was performed on one of two platforms: Agilent 244K oligonucleotide array, as described [START_REF] Stankiewicz | Genomic and Genic Deletions of the FOX Gene Cluster on 16q24.1 and Inactivating Mutations of FOXF1 Cause Alveolar Capillary Dysplasia and Other Malformations[END_REF] and the Affymetrix SNP 6 0 genotyping platform, as described [START_REF] Bochukova | Large, rare chromosomal deletions associated with severe early-onset obesity[END_REF] . An initial pilot study was performed using the Agilent 244K array; subsequently, the Affymetrix SNP 6 0 array was used in order to facilitate comparison with a cohort of 7,366 population controls, which had been assayed using the same platform. These control individuals were of European ancestry and were recruited from the WTCCC2 and GAIN (Genetic Association Information Network) consortia as described [START_REF] Bochukova | Large, rare chromosomal deletions associated with severe early-onset obesity[END_REF].

Affymetrix SNP 6 0 data were analyzed using Affymetrix powertools and Birdsuite software.

CNV calls made using Agilent software were converted to genotyping calls in order to enable comparison with the control cohort. We used a cut-off for calling CNVs of 200 Kb.

Validation of results

qPCR on patient DNA was used to confirm deletions or amplifications at specific genomic loci.

Results were normalized to ACTB DNA levels, which served as a control. Primer and probe sequences are shown in Table 4. Experiments were performed as described earlier [START_REF] Rad | CD25+/Foxp3+ T cells regulate gastric inflammation and Helicobacter pylori colonization in vivo[END_REF] on the ABI PRISM 7900HT Sequence Detection System (Applied Biosystems). Parental samples were not available for analysis; thus, information about whether a given copy number variant arose de novo was not available.

The Corresponding Author has the right to grant on behalf of all authors and does grant on behalf of all authors, an exclusive licence (or non exclusive for government employees) on a worldwide basis to the BMJ Publishing Group Ltd to permit this article (if accepted) to be published in JMG and any other BMJPGL products and sublicences such use and exploit all subsidiary rights, as set out in our licence (http://group.bmj.com/products/journals/instructions-for-authors/licenceforms). Affymetrix 6 0 data for three patients with deletions at 17q12 (Case 7 was not run on this platform). Log2 ratios for the three samples are highlighted in dark red, with the other samples from the same genotyping plate shown in black. The segmental duplication structure, taken Please refer to 8 for explanation of VCUAM classification of Mullerian abnormalities. Briefly, V5b = complete atresia of the vagina; C2b= bilaterlal aplasia of the cervix; U4b= aplasia of the uterus; A= adnex malformations with A1b: bilateral tubal malformation but normal ovaries; A2a: unilateral tubal hypoplasia; M= associated malformation, R: renal system; S: skeleton; C: cardiac; N: neurological; 0= normal).

existence of a previously undescribed 4 . 6

 46 Mb deletion on chromosome 2p. Clarification of the possible contribution of these two individual CNVs to abnormalities of Mullerian development must await further examples.

Figure 1 Copy

 1 Figure 1

Figure 2

 2 Figure 2

Table 2

 2 Summary of the clinical findings in patients with copy number variants of genomic disorder type and isolated or syndromic Mullerian aplasia

	Case	Locus of	Age	Classification	Findings	Isolated/	Renal	Cervical/vertebral	Cranio-	Cognitive	Growth	Other
		Deletion	(yrs)			syndromic			facial
	1	16p11.2	30	V5b C2b U4b A0	blind-ending	Syndromic	Normal	hypoplasia of the	No clinical	Moderate	Normal	Epilepsy,
				MSN Type II/ MURCS	vagina, rudimentary uterus, tubes and ovaries		(computer tomography and laparotomy)	wrist	indications	disturbed psychomotor development	(no exact measures)	moderate bilateral hearing loss
					normal
	2	16p11.2	20	V5b C2b U4b A0 M0	blind-ending	Isolated	Normal	Normal	No clinical	Normal	Ht (146 cm)
					vagina,		(computer		indications	psychomotor	Wt (47 kg)
				Type I	rudimentary		tomography			development
					uterus, tubes		and
					and ovaries		laparotomy)
					normal
	3	16p11.2	32	V5b C2b U4b A0 MR	blind-ending	Syndromic	Left atrophic	Scoliosis	No clinical	Normal	Normal
				Type II	vagina, rudimentary uterus, long		kidney		indications	psychomotor development	(no exact measures)
					uterus horns,
					tubes and
					ovaries
					normal
	4	16p11.2	18	V5b C2b U4b A0 M0	blind-ending	Isolated	Normal	Normal	No clinical	Normal	Ht (161 cm)
				Type I	vagina,		(laparascopy		indications	psychomotor	Wt (74 kg)
					rudimentary		and			development
					uterus, tubes		laparotomy)
					and ovaries
					normal
	5	17q12	24	V5b C2b U4b A0 M0	blind-ending	Isolated	Normal	Normal	No clinical	Normal	Ht (171 cm)
					vagina,		(computer		indications	psychomotor	Wt (56 kg)
				Type I	rudimentary		tomography			development
					uterus, tubes		and
					and ovaries		laparotomy)
					normal
	6	17q12	37	V5b C2b U4b A0	blind-ending	Syndromic	Left kidney	Pelvic	No clinical	Normal	Ht (173 cm)	Diabetes type 2,
				MRSC	vagina,		agenesis	misalignment	indications	psychomotor	Wt (70 kg)	cortico-
					rudimentary		with absent	resulting in		development		pulmonary
				Type II/ MURCS	uterus, tubes		ureter	different leg length				septal defect
					and ovaries
					normal

Table 4 .

 4 Primer and probe sequences used for quantitative PCR

	Region	Primer/Probe	Sequence
	1q21.1	Pex11b fw	CACGCTGATGTGCTTGTGATG
		Pex11b rev	TTTGACAATGATGAGGCCTGAA
		Pex11b probe	TGCGGCCTGCACTGGCCC
	1q21.1	CD160 fw	TTGAACCCAGGAGTCCACAAG
		CD160 rev	ACAGACGGCGGGAAACTCTT
		CD160 probe	CCAGGACCGCCTCCGAAGGTG
	2q11.2	NCAPH fw	GGGCGGCCTCCTCCTT
		NCAPH rev	TTCCAGGAAAACCACCATTTTAA
		NCAPH probe	CTCCTAAAGCGTGCTCGGTGTCTCTCC
	2q11.2	Kiaa1310 fw	GGTGCTGGCCAAGCAAGT
		Kiaa1310 rev	GCTCACCGACCCCAAGGT
		Kiaa1310 probe	TACTGTCTGTCCACGCGAGGTCTTTCTG
	16p11.2	MAZ fw	GGCTCAAAGGGCCCAATAG
		MAZ rev	CCTCCCTGTGCCCAGAAGT
		MAZ probe	AGGGATGCCCATGTACCACTCAGGC
	16p11.2	Tmem219 fw	GCACCCCACTTGGAAGCA
		Tmem219 rev	TGAGGCTCGCGGACTTTAA
		Tmem219 probe	TCAGATCTTGGCCCTACCCCTCCTGT
	17q12	LHX1 fw	CATGTGCCTGGGAAGAAAGG
		LHX1 rev	TGCCCTGTCTCTTCCAAGCT
		LHX1 probe	AGCCTGACTCGGCCCAGAAGCC
	17q12	Dusp14 fw	TCTGGTGCATGGATAGAAGCAA
		Dusp14 rev	CCACCGCAGAGAAAGACTCAA
		Dusp14 probe	TGACTTTCAGCGATGCCAAGTGTCCA
	2q11.2	Gnaz fw	GCCTGGTAGAGAGGTCTGTCTTG
		Gnaz rev	GGGAAATCACTTGGGCAGAA

*Probes were dual-labelled with FAM (5´) and TAMRA (3´)

Acknowledgements

The authors thank all of the patients who participated in this study and Depts of Obstetrics and Gynecology in Tübingen (Prof. D. Wallwiener), Heidelberg (Prof. T. Strowitzky) and Wuppertal (Prof. Hucke), Germany for help in recruitment of patients. C.S.-S., M.S., V.M., T.F., M.H., N.H., and N.P.C. are funded by the Wellcome Trust. C.S.-S. is the recipient of an Intermediate Clinical Fellowship from the Wellcome Trust. R.S. was funded by the Deutsche Forschungsgemeinschaft (DFG-STR923/2-2). The authors wish to thank Cordelia Langford and Sarah Widaa of the Microarray facility, and Sarah Edkins and Emma Gray of the DNA collections group at the Sanger Institute for assistance with sample processing and genotyping.

Conflicting Interests Statement

The authors declare no conflict of interest.

The Corresponding Author has the right to grant on behalf of all authors and does grant on behalf of all authors, an exclusive licence (or non exclusive for government employees) on a worldwide basis to the BMJ Publishing Group Ltd to permit this article (if accepted) to be published in JMG and any other BMJPGL products and sublicences such use and exploit all subsidiary rights, as set out in our licence (http://group.bmj.com/products/journals/instructions-for-authors/licenceforms). from the UCSC genome browser, is shown. Protein-coding genes are indicated by dark blue lines, with LHX1 and HNF1B highlighted in red.

Figure 3

Quantitative assessment of deletions and amplifications at indicated genomic loci.

Quantification was performed on DNA of 4-6 control persons (red bars) and affected individuals (blue bars). Primer and probe sequences were chosen to amplify intronic regions of the indicated genes. Data were normalized to results obtained at the ACTB locus and are presented as relative values.

Tables