The problem of efficiently compressing a large number, L, of N dimensional signal vectors is considered. The approach suggested here achieves efficiencies over current pre-processing and Karhunen-Loeve techniques when both L and N are large

Elhay, S

Gray, D

Howlett, P

Moran, William

Torokhti, A

Flinders Academic Commons

International Symposium on Signal Processing and its Applications, ISSPA, Gold Coast, Australia, 25-30 August, 1996. Organized by  the Signal Processing Research Centre, &UT, Brisbane, Australia. A NEW TECHNIQUE FOR MULTIDIMENSIONAL SIGNAL COMPRESSION A .  T o r ~ k h t i ~ ~ ~ * ~  D.Gray',+ S. Elhay4 P. Howlett3 W. Moran'r2 Cooperative Research Centre for Sensor Signal and Information Processing, The Levels, S A  5095, Australia Flinders University of South Australia, Adelaide, S A  5001, Australia University of South Australia, The Levels, SA 5095, Australia University of Adelaide, Adelaide, SA 5005, Australia 1. ABSTRACT The problem of efficiently compressing a large number, L ,  of fl dimensional signal vectors is considered. The approach suggested here achieves efficiencies over current preprocess- ing and Karhunen-LoCve techniques when both L and N are large. Preprocessing and partitioning techniques are first a p  plied to  the L x hf data matrix 3 to reduce the database to  a manageable number of subblocks of lower dimension. Within each subblock an iterative chain a p p r o h a t i o n  is proposed that effects a transform at  each stage of the it- erative scheme. A particularly appealing transform, using prolate spheroidal sequences, is suggested. To evaluate a reduced dimensionality approximation for the expansion coefficients, the approach used in the orthog- onal Procrustes problem solution is combined with an iter- ative interlacing technique due to  Daugavet for factorizing matrices. 2. INTRODUCTION Efficient compression of a reasonably large number, L, of signals, each represented by an N-vector is an important problem in the area of multidimensional signal processing. In the image processing case, each component of such a signal vector represents an image pixel. Typical applica- tions indude handwritten characters and face and object recognition. A novel application is the storage of ionogram back-scatter data  from over the horizon radar. W e  wiU deal with two dimensional signals, but the a p  proach we present easily extends to  the transmission and storage of signals with abitrary dimension. Let 3, the data  array formed by the vectors above, be the L x N matrix t o  be transmitted or stored. The choice of compression algorithm is determined, es- sentially, by the following competing criteria: the compres- sion ratio or the transmission cost, the quality of the data reconstruction (QDR), and the computational complexity. The methods commonly used in practice for data (image) compression or coding are based on computing a truncated Singular Value Decomposition (SVD) of the matrix 3 to approximately represent the data  array 3 as the product of two matrices, say, r and T, 3 = r x T, with smaller sizes L x M and M x N respectively, M < min(L,N). In the case of a large dimensionality of the matrix 7 ( t y p  ically L,  N = O(lO') - 0(106)) the processing required to  compute a truncated SVD can be prohibitive. Another difficulty is the data  transmission cost. Trans- mission of two matrices r and T is cheaper than the initial data array 3 transmission. However, for the case of the large L and N discussed above, the cost of transmitting matrices r, can be much greater than the cost of their evaluation, and may become prohibitively high. Thus, it is natural to strive to determine a technique which allows a further reduction in size of the transmitted or stored ma- trices while maintaining a reasonable QDR. The motivation of this paper comes from a need to reduce the overall cost of data array transmission. The method proposed allows a reduction in the total cost of the process of data compression, transmission/storage and reconstruction whilst balancing the criteria discussed earlier. To achieve this object we propose: a a new generic scheme of data array compression based on a combination of several preliminary stages and a final stage; m new discrete transforms for the preliminary stages; 0 a new compression algorithm for the final stage. 3. NEW DATA COMPRESSION METHOD Our approach to the problem solution is illustrated by the block-scheme shown in figure 1. I New genenc scheme I Figure 1. Block scheme of the signal compression method. The blocks (a), (b) and (c) in this figure mean the fol- (a): Preprocessing We first assume that the data  base 3 a has been preprocessed t o  remove redundanaes and similarities in the 3 by methods, for example, [6], and after that 0 has been partitioned into subblocks F ( 3 k ) ,  j = 1 : p ,  k = 1 : q, according to  schemes such as those consid- ered in [7]. After the above initial preprocessing stages the size of the subblock F ( J k )  F is 1 x n and their number is equal t o  P x q .  (b): Chain transform. The next step processes each of the pxq subblocks F ( J k ) .  A new generic scheme for carrying out this processing is proposed. The key idea of this generic scheme is a sequentid re- duction of dimensionalities of matrices in a matrix product which approximates each subblock F ( I k ) ;  these matrices are determined by a discrete transform which is applied a t  each lowing. of L signals (images): 162 International Symposium on Signal Processing and its Applications, ISSPA, Gold Coast, Australia, 25-30 August, 1996. Organized by the Signal Processing Research Centre, QUT, Brisbane, Australia. step of the reduction scheme. This part can be mathemat- ically represented as a solution of the problem with Tm,Tm-l ,..., TI given where F = F ( J k ) ,  W, - W k k )  E 72"q"' and T, = elk) E 7 L q 1 - 1 x q l , q o  = n,z = 1 : m, g m  < ... < ql < min (I,  n ) .  By using the equality F - WmT:. . . TT = ( F  - W1 TT)  + (WI - WZ T:)TT +... + (Wm-] - WmT:)T:-I ... TF (2)  (3)  problem (1) reduces to  1lWk-i - WkTTII -* mi: for k = 1 ,2 ,  ...: m and WO = P.We call this the iterative scheme of chain approximation. Under certain orthonormality conditions, the solution o b  tained is wk = Wk-lTk SO that Wm = FTiT2 ... Tm. Dia- grammatically this block (a) of the scheme above is repre- sented by a chain as follows The discrete transforms used in this process could be con- structed using prolate spheroidal sequences [l] or Cheby- shev polynomials. We briefly discuss the reason why discrete prolate spheroidal sequences are of interest here. They are defined as follows. Fix an integer N > 0 and some positive real number y E (0, i). We let T = T 7 , ~  be the linear operator on 1' given by where v(k) denotes the kth coordinate of v E 1'. This linear operator is symmetric with finite rank and so [I] it has N real eigenvalues Ay7, A;,. . . , > 0. Let v y  (k = 0,1,. . . , N - 1) be the corresponding eigenvec- N-l N-l - -  .. - tors normalised so that v y ( m ) '  = 1, v r P T ( m )  2 m=O mmO N-I 0, x ( N  - 1 - 2 m ) v y ( m )  3 0. These eigenvectors are called the discrete prolate spheroidal sequences. They have many interesting properties (loc. eit.), some of which will be of significance to  us. Their restrictions to the coordi- nates (0,1,2,. . . , N - 1) are the eigenvectors of the matrix whose (m,  n)th component is m=O sin 27ry(n - m) z(n - m) (n,m = 0 , 1 , 2  ,..., N - 1) An important property, the orthonormality of this index- limited collection of vectors, stems from the symmetry of the matrix. The nnlimited vectors are also orthogonal. Apart from these properties the sequences satisfy various symmetry laws and a difference equation. Of more interest to  us is that they (the unlimited sequences) are the most concentrated orthogonal sequences v E 1' whose Fourier series V ( f )  = C v ( n ) e a v n f  vanishes outside the interval n ( - 7 . 7 )  (modulo 1). That is, given the normalization above. the1 have the smallest 1' norm. The limited sequences are. on the other hand, those orthogonal sequences of length A- whose spectrum is most concentrated in the interval (-y. y )  in L 2 ( - i ,  +). The most interesting fact for our purposes concerns the set G, of sequences v in 1' which have lv(n)lZ _< E and where I = { 0 , 1 , 2 , .  . . , N - 1) . We look for the small- est M such that there is some M dimensional subspace D of 1' with the property that every member of G, can be approximated on the coordinate set I to  within 1) by a member d of D, that is, x ( g ( n )  - d(n))2 < 7. This M depends on N,  7, c and 1. For q > E > 0 the limit of M / N  as N - oc) is 27 and the discrete spheroidal se- quences vf*' (k = 0,1 ,2 , .  . . , [2Ny( l  - T)]) can be used as an orthogonal basis for an optimal space D in this sense. Thus for large N, if the sequences we are interested in are approximately of bandwidth 7 and are very small outside the index set I, as described above, then these sequences can be essentially described in terms of the first [27N] coef- ficients of their expansions in terms of the discrete prolate spheroidal sequences. In our context if the rows of the ma- trix which we are trying to  approximate have this property it makes sense to  use as the basis matrix one whose rows are the discrete prolate spheroidal functions. In effect the property is that there is relatively little high frequency (that is, close to *L) in the sequences. In the simplest case of using just one "link" in the trans- form chain above, i.e. m = 1 in (l), this approach was initially proposed in [2] and a particular case was consid- ered in [3]. (c): Interlacing iteration algorithm. For the find step of this generic scheme, a new matriz compression algorithm which combines the approach used in the solution of the orthogonal Procrustes problem [4] and the idea of interlac- ing iterations [5] is proposed and justified. This algorithm, which we term the Interlacing Iteration Algorithm, IIA, re- duces an evaluation of SVD for the matrix of interest to  eigen-decomposition evaluations for several special matri- ces with much lower dimensionalities as follows. For any matrices C = dd"' E 'RlXqm+l, and B = B$k) E 'Rqmxqm+l where qm+l < qm, we have n E I  F-CmBzTs ... TT = ( F - W ~ T ~ ) + ( W I - W Z T , ' ) T ~ + - - -  +(Wm-l- WmT;)T:,1. . . TT + (Wm - Cm B z ) T s .  . TT Then, taking into account (2), (3), the problem subject to  is reduced to  the problem sequence (3) considered above and the Karhunen-L&ve problem B ~ B  = I subject to  B ~ B  = I  (5) 163 International Symposium on Signal Processing and its Applications, ISSPA, Gold Coast, Australia, 25-30 August, 1996. Organized by  the Signal Processing Research Centre, Q UT,  Brisbane, Australia. where W = W" E 72' ",. For its solution, we first consider the related problem for a symmetric matrix S = WTW E 7LqmXqm: find X E 7Lqmxs and Y E 7Lqmxs with s << qm and r ( X )  = r ( Y )  << r ( W )  satisfying the problem 11s - XY' I I  - min (6) x ,  Y subject to Y T Y  = I ( 7 )  Solution of this problem is given by the following. position of the symmetric matriz A = A A  Lemma 1. Let A = X T S  and let the epenvalue decom- E 72''" be U B U T  = A ,  U T U  = U U T  = I .  Then X and Y satisfying the equations X = S Y  Y = SX(UQ-1 '2UT)  solve (6) and (7) .  Proof is based on using the Procrustes problem solution To evaluate X and Y which satisfy these equations we construct an iterative procedure as follows: Given X O  and S; for a = 0,1,2, ... put 141. A ;  = X:S (8) A, = AiAT (9) UiB,UiT = A ,  (10) x+1 = SX;(UiE);"2U,T) (11) x,+1 = sx+1 (12) We point out that the fact the eigen-decomposition of the small s x s, symmetric matrix A, is easy to  compute is central to the appeal of IIA. Choice of the initial itemtion X O  requires special atten- tion. Suppose ranks = r so that s 5 r 5 qm. Put  X O  = Sa E Rqmx' where S, consists of s linearly inde- pendent columns of S. Let the SVD of S, be sa = @.INT with ip E Rqmxr,C E a'., it? E R'". Then xo = QP (13) where P E RrXa. If the first s rows of the matrix P in (13) are linearly independent then the corresponding matrix xo is denoted- X;. Lemma 2. Let the pair X ,  Y be a solution of the problem (61.f7). Let 1 > am and let S = W T W  have r a n k s  = r .  (If  . , . I  I qm > I we define S = W W T ) .  Denote b y  XI, XZ, ...,A, the eigenvalues of S which are such that 1x11 2 1x21 1 ... 2 I X r I  > 0 Suppose that s < z,Xo = X ;  and that As+] # A,. Then, the following limits m 's t  x = lim X.+I ,  Y = lim x+1 1-00 ,-OD and we have IIX - X.+III < cI(X.+I/X')z('+I) IIY - X+III < cz(X,+l/x,)z('+l) with constants c1. cz which are independent of A,. X,+1 and 1 .  Proof is derived from [ 5 ] .  The lemmas 1,2 allow us to state the following theorem which includes solutions of the problem (4), (5)  and cor re  spondmg error estimations. Theorem 1. Let Z = X T Y  where X ,  Y is a solution of (6),(7) and let the eigen-decomposition of S be s =  GAG^, G ~ G  = G G ~  = I  with r a n k s  = r andA = A, E 72"'. L e t s  < r ,  and let Z = VRVT, V E R S x "  orthogonal, be the eigen-decomposition of the symmetric matriz Z .  Then, R = A ,  and the solution required is defined by the relations C =  W B ,  B = Y V  Moreover, CTC = Theorem 2. Let the conditions of Theorem 1 be ful- filled and let K+i be defined b y  the eigen-decomposition V,+I D,+I V,T1 = Z,+I where Z,+I = X z l  X+I.  If A,+] # A,, then for  C.+I = WB,+i. (14) B,-LI = Y,+iV,+1 (15) we haoe IIC - C,+IlI < c3(x,+I/Xs)z('+') I IB  - ~ i + l l l <  C~(X,+I/X,)~('+') with constants c3, c4 which are independent of A,, A,+] and Proofs of these theorems can be established drawing on the results [SI. On the basis of the above, the interlacing iteration a1g.o- rithm, IIA, for solution of the problem (4), (5) is as follows: 1. given W,S and X O ;  compute Y;+I and Xi+]  for i = 0 ,1 ,2 ,  .__ in accordance with (8)-(12); 2. when X,+I and Y,+I satisfy the stopping criteria IIXt+l - Xi11 5 e', IIX+I - Y,II 5 e'' with E',"' given, compute the small (s x s) dimension symmetric matrix Z;  = X:X for this i and the eigen-decomposition 1. zi = v; ni KT (16) where ~ ? i  = diag(X1, XZ, ..., A,} as in Theorem 1. 3. Finally, compute the matrices Bi = X K ,  Ci = W B i  (17) Note that only C, and B;, of much reduced dimension- alities, are transmitted and this leads t o  a correspond- ing reduction of the transmission costs. Let us consider the computational complezity of IIA. S u p  pose, l > qm in (4) ,  (5). It follows from (8)-(12) and (14), (15), that each iteration of IIA takes O(sqm) flops to compute the matrices Ai, X+1 and Xi+]. It also requires O(iq,q,+l) flops to evaluate the matrix Ci. These counts determine the dominant part of the computational com- plexity of IIA. By comparison, using a QR-type algorithm to solve the problem (4), (5) for C and B requires 0(413/3)  flops per iteration. Also note that IIA has the high degree of convergence. 164 International Symposium on Signal Processing and its Applications, ISSPA, Gold Coast, Australia, 25-30 August, 1996. Organzzed by  the StgnaI Processing Research Centre, Q UT, Brisbane, Australia. 0.8 4. NUMERICAL EXAMPLES To demonstrate the viability of the IIA method we took some electro encephalographic (EEG) data in which there are 512 rows, each with 167 readings taken from one EEG sensor. This data was normalized to have maximum mod- ulus unity. The rows are grouped in sets of 32, each group representing one test with 32 simultaneously sampled sen- sors. We applied the compression technique above to the first 3 = 32 rows, the first j = 128 rows and then to all j = 512 rows. For each choice of row dimension we ran IIA with s = 4,12,16 to give a range of compression ratios. On each test, we took the first s columns of the matrix S as a starting matrix XO and ran the iteration (8)-(12) until the relative difference, between successive approximations, of Frobenius norm 11s - X Y T l l ~  wbs smaller than All the calculations were done using IEEE double precision arithmetic with a machine epsilon of about 2 x and the Frobenius norm is used throughout. The results of these illustrative examples (not meant to be a thorough numerical investigation) are summarized in the Table 1 below. Column 1 of this table shows the number of rows in the sample and Column 2, the number s. The next two columns, 3 and 4, show the ratios - 1 respectively, where p is the number of iterations required for convergence and W A ~  is the approximation to W given by IIA. Column 5 shows the compression, C p  = (< C > + < B >)/ c W >,' where < . > denotes the number of elements in the argument matrix, and Column 6 shows the number p of iterations for convergence. We observed that, universally, the major reduction in the norms (18) occurs in the first iteration: subsequent itera- tions only produce a small decrease in the accuracy of the approximation. This is illustrated in Table 2, where the error a t  each of the four iterations needed for convergence, is shown. The compression for this example has reduced the size of the approximation to 16.5% of the original size. Clearly, the improvement obtained depends on the choice of initial matrices XO and YO. However, even for arbitrar- ily chosen XO and YO we saw very similar convergence be- haviour. In Figure 2, we show the corresponding EEG trace and its approximation (i.e. the reconstruction after com- pression by IIA) on the same axes. Table 1. Numerical examples: error ratios, compression and number of iterations Table 2. 128 Error reduction for rows with s=12 o 20 40 so 80 io0 120 140 160 teo -0.4 167datapoiw Figure 2. A typical EEG trace and its reconstruction 5. ADVANTAGES AND SIGNIFICANCE In comparison with the known methods 11-51 - -  The specific generic scheme of the signal compression proposed makes it possible to  considerably cut the transmission cost; The new discrete transforms suggested lead to a high quality of the data compresed reconsmction; The new approach of uncorrelated compression based on a specific combination of the ideas [4], p. 582. and [SI has a high speed of convergence and consequently allows us to  economize the computational expenditures needed to  reach the desired accuracy. REFERENCES [I] D. Slepian, "Prolate spheroidal wave functions, Fourier analysis, and uncertainty - V: the discrete case", The Bell System Technical Journal, vol. 57, no. 5, pp. 1371-1430, 1978 [2] A. Torokhti, "On polynomial synthesis of nonlinear systems", Electron. Modelling, .I. of USSR Acad. of Sc., vol. 11, no. 4, pp. 28-34, 1989. In Russian. [3] H. Murase, M. Lindenbaum, "Partial eigenvalue de- composition of large images using spatial temporal adaptive method", IEEE Tmns. Image Processing, vol. 4, no.5, pp. 620-629, 1995. [4] G.H. Golub and C.F. Van Loan, M a t h  computa- tions, Johns Hopkins University Press, 1989. [5] V.A. Daugavet and P.V. Yakovlev, "Mean-square' approximation of a rectangular matrix by matrices of lower rank", Zh. Vychisl. Mat. a' Mat. Fiz., vol. 29, no. 10, pp. 1466-1479, 1989. In Russian. [SI P.J.E.Vermeulen and D.P.Casasent, "Karhnnen- L o k  techniques for optimal processing of time- sequential imagery", Optical Engineering, 30 (1991), 415-421. [7] J. Vaisey, A. Gersho, "Image compression with vari- able block size segmentation", IEEE Tmm. on Sag- nal Processing, 40, (1992), 2040-2060. 165 

A New Technique for Multidimensional Signal Compression

http://dspace.flinders.edu.au/xmlui/bitstream/handle/2328/25773/Torokhti%20New%20Technique.pdf?sequence=1

A New Technique for Multidimensional Signal Compression

Abstract

Similar works

Full text

Available Versions

Flinders Academic Commons