CA2851028A1 - Method and system for compressed sensing of physiological signals - Google Patents

Abstract

A method and system are provided for performing compressed sensing seismocardiography. An input vector x is received and it is encoded with a sensing matrix .PHI.. The encoded data y is decoded by finding a sparse solution to argmin, s. t. y = .PHI..NU..
v¦¦.NU.¦¦O
An approximation x of the input vector y is reconstituted by calculating x =
.PHI..PSI., where .PSI. is a fixed-basis transformation matrix that is selected based on the type of signal that is being measured. Particular bases for use with seismocardiography are provided.

Description

METHOD AND SYSTEM FOR COMPRESSED SENSING OF PHYSIOLOGICAL
SIGNALS
Technical Field [0001] This disclosure relates to signal sampling, and in particular to methods and systems for compressed sensing.
Background

[0002] Some signals-processing systems sample information-dense signals, often requiring that large quantities of information be acquired and/or processed in the course of sampling. This may be costly, infeasible, or otherwise undesirable, particularly when the available bandwidth for acquiring, transmitting or processing samples is limited. As a consequence, there is a need for ways to represent signals using smaller amounts of data, while preserving important characteristics of sampled signals without unacceptably degrading the accuracy of the system's measurements of the sampled signal.

[0003] It is well known that the sampling rate of a signals-processing system can be reduced to the Nyquist rate (i.e. twice the maximum frequency of the signal being sampled) without a loss in accuracy. However, even at the Nyquist rate, some systems still require the acquisition and/or processing of substantial quantities of data.
Examples of such systems in the medical field include magnetic resonance imaging (MRI) systems, ultrasound imaging systems, ballistocardiogram (BCG) systems, and seismocardiogram (SCG) systems.

[0004] One approach to reducing the sampling rate below the Nyquist rate (and thereby reducing the amount of data that must be acquired to reconstruct a signal) is compressed sensing. Compressed sensing allows signals with particular characteristics ¨
namely, signals that are sparse ¨ to be sampled at rates below the Nyquist rate without losing any information about those signals. In practice, perfectly reconstructing signals sampled with compressed sensing can be difficult (i.e. NP-hard), but good approximations of the original signal can be reconstructed using numerical analysis.

[0005] Compressed sensing takes advantage of signal sparsity to reduce the number of measurements required to reconstruct a signal. Sparsity can be defined as follows: The sparsity of an Nxl vector p is given by its number of non-zero elements. This quantity is also known as the Lo norm, and can be denoted as II 110. A vector with low sparsity is said to be sparse. The more sparse that a signal is, the more effective compressed sensing tends to be at compressing the signal.

[0006] One of the challenges in performing compressed sensing for any particular application is obtaining a sparse representation of the relevant signal. Some signals are naturally sparse. For example, electrocardiogram (ECG) signals tend to be relatively sparse in the time domain; an example of such a signal is shown in Figure 4.
However, other signals, such as SCG signals (an example of which is also shown in Figure 4), tend to be significantly less sparse in the time domain, reducing the efficacy of compressed sensing when applied directly to those signals.

[0007] Some signals that are not naturally sparse may be transformed into equivalent signals that are sparse, thereby improving the efficacy of compressed sensing.
For example, truncation has been applied to BCG signals so that values below a threshold are rounded down to zero. Truncation reduces the accuracy of the measurements, resulting in lossy compression. Whether a particular transformation is likely to induce sparsity of a signal, and whether that transformation will have a detrimental effect on accuracy, is often highly dependent on the characteristics of the signal.
Summary

[0008] The present disclosure provides several aspects and embodiments, including the following:

[0009] Aspect 1 provides a method of performing compressed sensing seismocardiography. The method comprises receiving an input set of data representing a seismocardiogram signal, obtaining a sensing matrix, determining an encoded set of data based on a product of at least a portion of the input set of data and the sensing matrix, obtaining a transformation matrix corresponding to the sensing matrix, and obtaining a solution matrix based on a product of the sensing matrix and the transformation matrix.
The solution matrix satisfies the restricted isometry property. The method also comprises determining a sparse interim set of data so that a product of the sparse interim set of data and the solution matrix corresponds to the encoded set of data and determining an output set of data based on a product of the transformation matrix and the sparse solution so that the output set of data approximates the input set of data.

[0010] Aspect 2 provides a method according to aspect 1 wherein the sensing matrix has a size M xN where N is a size of the input set of data and M is less than N.

[0011] Aspect 3 provides a method according to any one of aspects 1 and 2 wherein the sensing matrix is constructed from independent and identically distributed Gaussian distributions.

[0012] Aspect 4 provides a method according to any one of aspects 1 to 3 wherein obtaining the sensing matrix comprises retrieving the sensing matrix from a memory.

[0013] Aspect 5 provides a method according to any one of aspects 1 to 3 wherein obtaining the sensing matrix comprises generating the sensing matrix.

[0014] Aspect 6 provides a method according to any one of aspects 1 to 3 wherein obtaining the sensing matrix comprises receiving the sensing matrix from a user.

[0015] Aspect 7 provides a method according to any one of aspects 1 to 6 wherein the transformation matrix has a size NxN where N is a size of the input set of data.

[0016] Aspect 8 provides a method according to any one of aspects 1 to 7 wherein the transformation matrix and the sensing matrix are incoherent.

[0017] Aspect 9 provides a method according to any one of aspects 1 to 8 wherein the transformation matrix has a fixed basis.

[0018] Aspect 10 provides a method according to aspect 9 wherein the fixed basis comprises a PDCT so that the kth column of the fixed basis is based on the construction given by:

[1, COS (¨kn.) COS ((N 1)kir T
N N
where N is a size of the input set of data.

[0019] Aspect 11 provides a method according to aspect 10 wherein the fixed basis has normalized zero-mean columns.

[0020] Aspect 12 provides a method according to aspect 9 wherein the fixed basis comprises the Haar wavelet transform.

[0021] Aspect 13 provides a method according to aspect 9 wherein the fixed basis comprises the Daubechies 2 wavelet transform.

[0022] Aspect 14 provides a method according to aspect 9 wherein the fixed basis comprises the Daubechies 5 wavelet transform.

[0023] Aspect 15 provides a method according to aspect 9 wherein the fixed basis comprises the Symlet 2 wavelet transform.

[0024] Aspect 16 provides a method according to aspect 9 wherein the fixed basis comprises the Symlet 4 wavelet transform.

[0025] Aspect 17 provides a method according to aspect 9 wherein the fixed basis comprises the Bioorthogonal 4.4 wavelet transform.

[0026] Aspect 18 provides a method according to aspect 9 wherein the fixed basis comprises the Bioorthogonal 5.5 wavelet transform.

[0027] Aspect 19 provides a method according to aspect 9 wherein the fixed basis comprises the Reverse Bioorthogonal 4.4 wavelet transform.

[0028] Aspect 20 provides a method according to aspect 9 wherein the fixed basis comprises the Reverse Bioorthogonal 5.5 wavelet transform.

[0029] Aspect 21 provides a method according to any one of aspects 1 to 20 wherein obtaining the transformation matrix comprises retrieving the transformation matrix from a memory.

[0030] Aspect 22 provides a method according to any one of aspects 1 to 20 wherein obtaining the transformation matrix comprises generating the transformation matrix.

[0031] Aspect 23 provides a method according to any one of aspects 1 to 20 wherein obtaining the transformation matrix comprises receiving the transformation matrix from a user.

[0032] Aspect 24 provides a method according to any one of aspects 1 to 23 wherein determining a sparse interim set of data comprises finding a sparse solution to:
argminilvilo , s.t. y =
where y is the encoded set of data and 4' is the sensing matrix.

[0033] Aspect 25 provides a method according to aspect 21 wherein finding a sparse solution comprises performing basis pursuit.

[0034] Aspect 26 provides a method according to aspect 21 wherein finding a sparse solution comprises performing matching pursuit.

[0035] Aspect 27 provides a method according to aspect 21 wherein finding a sparse solution comprises performing orthogonal matching pursuit.

[0036] Aspect 28 provides a method according to any one of aspects 1 to 27 wherein the seismocardiogram signal is a tri-axial seismocardiogram signal.

[0037] Aspect 29 provides a method according to any one of aspects 1 to 28 wherein the input set of data is segmented according to a windowing scheme.

[0038] Aspect 30 provides a method according to any one of aspects 1 to 29 comprising obtaining a second sensing matrix, using the second sensing matrix in place of the sensing matrix, and obtaining a second solution matrix based on the second sensing matrix. The second solution matrix satisfying the restricted isometry property. The method also comprises using the second solution matrix in place of the solution matrix.

[0039] Aspect 31 provides a method according to aspect 30 comprising determining an accuracy measure of the output set of data, wherein the second sensing matrix is obtained in response to the determination of the accuracy measure.

[0040] Aspect 32 provides a method according to any one of aspects 30 and 31 wherein the second sensing matrix is obtained in response to user input.

[0041] Aspect 33 provides a method according to any one of aspects 30 to 32 wherein the second sensing matrix is obtained based on a desired compression ratio.

[0042] Aspect 34 provides a system for compressed sensing seismocardiography.
The system comprises an encoder for encoding an input set of data. The encoder is configured to determine an encoded set of data based on a product of at least a portion of the input set of data and a sensing matrix. The system also comprises a decoder for decoding the encoded set of data. The decoder is configured to determine a sparse interim set of data so that a product of the sparse interim set of data and a solution matrix corresponds to the encoded set of data and to determine an output set of data based on a product of a transformation matrix and the sparse solution so that the output set of data approximates the input set of data. The transformation matrix corresponds to the sensing matrix, the solution matrix is based on a product of the sensing matrix and the transformation matrix, and the solution matrix satisfies the restricted isometry property.

[0043] Aspect 35 provides a system according to aspect 34 wherein the sensing matrix has a size M xN where N is a size of the input set of data and M is less than N.

[0044] Aspect 36 provides a system according to any one of aspects 31 and 35 wherein the sensing matrix is constructed from independent and identically distributed Gaussian distributions.

[0045] Aspect 37 provides a system according to any one of aspects 31 to 36 wherein the transformation matrix has a size N xN where N is a size of the input set of data.

[0046] Aspect 38 provides a system according to any one of aspects 31 to 37 wherein the transformation matrix and the sensing matrix are incoherent.

[0047] Aspect 39 provides a system according to any one of aspects 31 to 38 wherein the transformation matrix has a fixed basis.

[0048] Aspect 40 provides a system according to aspect 39 wherein the fixed basis comprises a PDCT so that the kth column of the fixed basis is based on the construction given by:
kir (N ¨ 1)kg T
[1, COS , COS ( (-N) N )1 where N is a size of the input set of data.

[0049] Aspect 41 provides a system according to aspect 40 wherein the fixed basis has normalized zero-mean columns.

[0050] Aspect 42 provides a system according to aspect 39 wherein the fixed basis comprises the Haar wavelet transform.

[0051] Aspect 43 provides a system according to aspect 39 wherein the fixed basis comprises the Daubechies 2 wavelet transform.

[0052] Aspect 44 provides a system according to aspect 39 wherein the fixed basis comprises the Daubechies 5 wavelet transform.

[0053] Aspect 45 provides a system according to aspect 39 wherein the fixed basis comprises the Symlet 2 wavelet transform.

[0054] Aspect 46 provides a system according to aspect 39 wherein the fixed basis comprises the Symlet 4 wavelet transform.

[0055] Aspect 47 provides a system according to aspect 39 wherein the fixed basis comprises the Bioorthogonal 4.4 wavelet transform.

[0056] Aspect 48 provides a system according to aspect 39 wherein the fixed basis comprises the Bioorthogonal 5.5 wavelet transform.

[0057] Aspect 49 provides a system according to aspect 39 wherein the fixed basis comprises the Reverse Bioorthogonal 4.4 wavelet transform.

[0058] Aspect 50 provides a system according to aspect 39 wherein the fixed basis comprises the Reverse Bioorthogonal 5.5 wavelet transform.

[0059] Aspect 51 provides a system according to any one of aspects 31 to 50 wherein to determine a sparse interim set of data comprises finding a sparse solution to:
argminfivIlo s=t= Y =
where y is the encoded set of data and 4 is the sensing matrix.

[0060] Aspect 52 provides a system according to aspect 51 wherein finding a sparse solution comprises performing basis pursuit.

[0061] Aspect 53 provides a system according to aspect 51 wherein finding a sparse solution comprises performing matching pursuit.

[0062] Aspect 54 provides a system according to aspect 51 wherein finding a sparse solution comprises performing orthogonal matching pursuit.

[0063] Aspect 55 provides a system according to any one of aspects 31 to 54 wherein the seismocardiogram signal is a tri-axial seismocardiogram signal.

[0064] Aspect 56 provides a system according to any one of aspects 31 to 55 wherein the input set of data is segmented according to a windowing scheme.

[0065] Aspect 57 provides a system according to any one of aspects 31 to 56, wherein the system is configured to obtain a second sensing matrix, use the second sensing matrix in place of the sensing matrix, and obtain a second solution matrix based on the second sensing matrix. The second solution matrix satisfying the restricted isometry property. The system is also configured to use the second solution matrix in place of the solution matrix.

[0066] Aspect 58 provides a system according to aspect 57 wherein the system is configured to determine an accuracy measure of the output set of data and to obtain the second sensing matrix in response to the determination of the accuracy measure.

[0067] Aspect 59 provides a system according to any one of aspects 57 and 58 wherein the system is configured to obtain the second sensing matrix in response to user input.

[0068] Aspect 60 provides a system according to any one of aspects 57 to 59 wherein the system is configured to obtain the second sensing matrix based on a desired compression ratio.

[0069] Some embodiments provide one or more apparatus having any new and inventive feature, combination of features, or sub-combination of features as described herein.

[0070] Some embodiments provide a method or methods having any new and inventive steps, acts, combination of steps and/or acts or sub-combination of steps and/or acts as described herein.

[0071] Further aspects and example embodiments are illustrated in the accompanying drawings and/or described in the following description.
Brief Description of the Drawings

[0072] The accompanying drawings illustrate non-limiting example embodiments of the invention.

[0073] Figure 1 is a block diagram showing an example sensing system.

[0074] Figure 2 is a flow diagram of an example method for compressed sensing.

[0075] Figure 3 is a flow diagram of an example method for compressed sensing.

[0076] Figure 4 is a block diagram of an example seismocardiography (SCG) system.

[0077] Figure 4a is a plot of an example x-SCG signal from an example tti-axial SCG
signal.

[0078] Figure 4b is a plot of an example y-SCG signal from the example tri-axial SCG
signal of Figure 4a.

[0079] Figure 4c is a plot of an example z-SCG signal from the example tri-axial SCG

signal of Figure 4a.
100801 Figure 4d is a plot of an example ECG signal corresponding to the example tri-axial SCG of Figures 4a-4c.
[0081] Figure 5a is a set of plots of reconstructed signals from an example experimental application of an embodiment of the invention in which the compression ratio is set to 1.
[0082] Figure 5b is a set of plots of reconstructed signals from an example experimental application of an embodiment of the invention in which the compression ratio is set to 1/6.
100831 Figure 6a is a set of plots of results of an example experimental application of an embodiment of the invention in which a variety of compression ratios are used, measured using peak signal-to-noise ratio (PSNR).
[0084] Figure 6b is a set of plots of results of an example experimental application of an embodiment of the invention in which a variety of compression ratios are used, measured using average normalized root-mean-square deviation (nRMSD).
[0085] Figure 6c is a set of plots of results of an example experimental application of an embodiment of the invention in which a variety of compression ratios are used, measured using maximum nRMSD.
[0086] Figure 6d is a set of plots of results of an example experimental application of an embodiment of the invention in which a variety of compression ratios are used, measured using minimum tiRMSD.
Description 100871 Throughout the following description specific details are set forth in order to provide a more thorough understanding to persons skilled in the art. However, well known elements may not have been shown or described in detail to avoid unnecessarily obscuring the disclosure. The following description of examples of the technology is not intended to be exhaustive or to limit the system to the precise forms of any example embodiment.
Accordingly, the description and drawings are to be regarded in an illustrative, rather than a restrictive, sense.
Sensing System 100881 Figure lshows an example compressed sensing system 1. Compressed sensing system 1 comprises a sensor 10 that takes samples of a signal and provides raw data s to a signal conditioner 15. Signal conditioner 15 performs preprocessing on data s and provides the processed data t to analog/digital converter 20. Analog digital converter 20 represents processed data tin a suitable (digital) format for later processing by other components of compressed sensing system 1, such as encoder 30 and decoder 50. Analog/digital converter 20 provides input data x, derived from processed data t, to encoder 30.
Analog/digital converter 20 may provide input data x to encoder 30 directly or indirectly;
for example, additional preprocessing may be applied to input data x, such as windowing of input data x in the time domain and/or in other domains. In some embodiments, analog/digital converter 20 also (or alternatively) provides input data x to detector 25, which is discussed further below.
100891 Encoder 30 takes data x as input and, based on that input, provides encoded data y.
Encoded data y can be used to reconstruct input data x (or an approximation thereof), and is smaller than input data x at least some of the time. Thus, by transmitting and/or storing encoded data y in the place of data x, bandwidth and/or data storage can be saved. Decoder 50 takes encoded data y as input and provides output 2, which is equal to or approximates data x.
100901 In some embodiments, sensor 10 and encoder 30 are provided in a single device.
For example, a housing may contain both sensor 10 and encoder 30, which may be in communication so that readings from sensor 10 are provided to encoder 30.
Encoded data produced by encoder 30 may be transmitted from the housing to decoder 50, which may be provided in a remote device, and/or the encoded data may be stored locally on the device and later decoded by a local or remote decoder 50. As another example, sensor 10 may transmit its readings or a derivation thereof to encoder 30, which may be provided in a remote device where encoder 30 may encode the received information. Encoder 30 may then transmit the encoded information to a remote decoder 50 and/or store it locally for later decoding by a local or remote decoder 50.
[0091] In some embodiments, encoder 30 performs a linear map:
y = Ox, (1) where x is an Nxl vector representing a block of the input data, cf) is an MxN
sensing matrix, and y is an Mxl vector representing the output. In at least some cases, M < N, resulting in a compression ratio of MN <1.
[0092] Decoder 50 reconstructs x (or an approximation thereof) from y. In general, finding such a reconstruction is a problem with no known solution, given arbitrary y.
However, under certain conditions, a sparse solution can be found that is the same as x. In some embodiments, decoder 50 finds the sparse solution:
= argrninllvll0, s.t. y = (fiv, (2) where y and are as defined for Eq. (1). Otherwise put, the sparse solution 2 is the sparsest vector that can be found that satisfies the equation y = 02. Once a sparse solution has been found, decoder 50 provides it as output.
[0093] In general, the sparse solution I is not guaranteed to be identical to the original data x. However, according to compressed sensing theory, 2 = x is guaranteed if sensing matrix q satisfies the restricted isometry property.
[0094] One of the known consequences of the restricted isometry property is, that for cases where satisfies the restricted isometry property, the size M of vector y has a lower bound:
M > K log¨K (3) where K is the sparsity of x and Nis the size of x. This implies that better compression ratios are possible for sparser input vectors x than are possible for less-sparse input vectors x.

[0095] Solving Eq. (2) optimally is an NP-hard problem. As a consequence, decoder 50 may find a sub-optimal solution for the sake of practicality. Several approaches for finding sub-optimal solutions of Eq. (2) exist, including basis pursuit, matching pursuit, derivatives thereof, and other approaches. This allows decoder 50 to find a sparse solution that approximates x to a reasonable degree of accuracy. 2 can then be subject to further processing and/or analysis.
Encoding Method [0096] Compressed sensing system 1 may perform, for example, compressed sensing method 100, as shown in Figure 2. In some embodiments, encoder 30 may perform encoding method 110 and decoder 50 may perform decoding method 150.
[0097] At block 115, encoder 30 receives data x. Data x is provided directly or indirectly to encoder 30 by sensor 10.
[0098] At block 120, compressed sensing system 1 obtains sensing matrix 0.
Sensing matrix .0 may be generated by encoder 30 after receiving data x, as shown in Figure 2, or before receiving data x. Alternatively, or in addition, sensing matrix 0 may be generated by compressed sensing system 1 in initialization method 205, as shown in Figure 3.
Sensing matrix 0 may, for example, be generated before or during signal acquisition. In some embodiments, compressed sensing system 1 is pre-programmed with one or more sensing matrices 0 for use in encoding method 110 and decoding method 150.
Compressed sensing system 1 may also, or alternatively, receive one or more sensing matrices 0 from a user.
[0099] Compressed sensing system 1 may use one sensing matrix .0 for one, some, or all vectors of data x received from sensor 10. Compressed sensing system 1 may use one sensing matrix 0. at all times or may use different sensing matrices .0 at different times, such as periodically, after processing each vector of data x, randomly, and/or at other times. Compressed sensing system 1 may switch to using a different sending matrix 0, for example, in order to achieve a different compression ratio. Such a switch may, for example, be in response to user instruction, a change in signal morphology, a change in the transformation matrix IP used by decoder 50, an increase or decrease in noise in the signal, and/or for other reasons.
[0100] Notably, sensing matrix 4) does not necessarily need to satisfy the restricted isometric property. Since x is not necessarily sparse, compressed sensing system 1 may induce sparsity on x. In general, we say that x is sparse under a transformation 111 if, for an NxN matrix W, there is an Nxl vector s such that x= T s, (4) where s is sparse (i.e., the sparsity K of s is less than the size N of x). If such a transformation matrix T can be found, then it is sufficient for sensing matrix 4) to be selected so that the product e= OW satisfies the restricted isometric property. This is due to the fact that Eq. (1) can be rewritten as follows:
y = (fix = 4)Ts = ()s, (5) allowing the decoder to solve for:
= argminlivilo , s.t. y = ev (6) instead of finding a solution to Eq. (2) directly.
[0101] In some embodiments, 4) and W are incoherent. Incoherence can be defined in several ways; the most common definition is column-wise incoherence, meaning that the magnitude of the inner product of columns 4)k and 4)i from matrices 4) and IP, respectively, is small for each pair (4)k, 4)i). Alternative definitions of incoherence may also be used (e.g. pair-wise, rank-wise, K-wise, and so on). Unless otherwise indicated, the following examples use the column-wise definition of incoherence.
[0102] In some embodiments, sensing matrix 4) is constructed randomly (or, more commonly, pseudo-randomly) from independent and identically distributed (i.i.d.) distributions which may be Gaussian distributions. Such matrices have very high probabilities of being incoherent with fixed-basis matrices. In embodiments where compressed sensing system 1 uses a fixed-basis transformation matrix 111, the incoherence between 0 and W means that e is likely to satisfy the restricted isometry property.
Sensing matrix 0 may be constructed using other methods; for example, in some embodiments sensing matrix 0 is constructed from a Gaussian circular distribution (which is not necessarily i.i.d.). Constructing sensing matrix 0 from a Gaussian circular distribution may, in some circumstances, result in increased signal degradation and/or reduced compression ratios during compression/decompression. Compression ratios are discussed at greater length below.
[0103] At block 125, encoder 30 encodes x by generating y according to Eq.
(1). At block 130, encoder 30 outputs y. The output of encoder 30 may be stored and/or transmitted, as described above.
Decoding Method [0104] At some point, encoded data y is received by decoder 50 at block 155.
Encoded data y is provided directly or indirectly to decoder 50 by encoder 30.
[0105] At block 160, compressed sensing system 1 obtains transformation matrix W.
Transformation matrix W may be generated by decoder 50 after receiving encoded data y, as shown in Figure 2, or before receiving encoded data y. Alternatively, or in addition, transformation matrix may be generated in initialization method 205, as shown in Figure 3. Transformation matrix W may, for example, be generated before or during decoding, and may be generated by compressed sensing system 1 or by another system and/or user. In some embodiments, compressed sensing system 1 is pre-programmed with one or more transformation matrices W for use in decoding method 150. However or whenever transformation matrix 11' is generated, it must be provided to decoder 50, e.g. at block 160.
[0106] The selection of the transformation matrix IP depends on the properties of the signal being encoded and decoded. That is, the transformation matrix W should be selected so that input vectors x describing the measure signal are likely to be sparse under W.

Different classes of signals may be better suited to certain transformation matrices W than other signals, and vice-versa. Specific methods of finding transformation matrices for particular classes of signals are discussed below.
101071 Compressed sensing system 1 may use one transformation matrix IP for one, some, or all vectors of encoded data y received from encoder 30. Compressed sensing system 1 may use one transformation matrix IP at all times or may use different transformation matrices LP at different times, such as periodically, after processing each vector of encoded data y, randomly, and/or at other times. For example, a different 'P may be used if the morphology of the signal changes, and/or if the characteristics of the input vectors change.
For instance, an SCG signal may have a particular morphology when measuring the activity associated with a particular range of heart rates, but that morphology may change if the heart rate dramatically increases or decreases, or if the physiology of the heart becomes distressed, or in other circumstances. If such a change is detected, compressed sensing system 1 may switch to a different transformation matrix IF that is known to be better-suited to the morphology of the new signal. As another example, encoder 30 may switch to a different sensing matrix , in response to which decoder 50 may switch to a different transformation matrix W.
101081 In some embodiments, compressed sensing system 1 selects a transformation matrix IP to use in encoding/decoding by detecting the signal before encoding the signal and determining which of a plurality of available transformation matrices µP
most effectively induces sparsity on that signal. Such detection and selection may be performed, for example, by a detector 25 that receives input data x from analog/digital converter 20, analyzes the input data x to determine which transformation matrix µP is most suitable for input data x (e.g. by applying each transformation matrix to input data x and determining the sparsity resulting from each application), selects a transformation matrix W accordingly, and passes input data x to encoder 30. In some embodiments, input data x is provided to both detector 25 and encoder 30 and they perform their functions in parallel;
the selection of transformation matrix IP is then communicated to decoder 50 prior to decoding.

[0109] In some embodiments, a transformation matrix IP may be selected from a given set of candidates (e.g. a family of wavelet bases known to be suitable for a particular class of signals). The candidates may, for example, be ranked, with one or more of the highest-ranked candidates being selected. Quantitative metrics may be used to rank the candidates;
for example, the candidates may be ranked according to measurements of the peak signal-to-noise ratio (PSNR) and/or normalized root-mean-square deviation (nRMSD) of decoded signals transformed with the candidate transformations. PSNR and nRMSD
are defined in greater detail below. In some embodiments, the basis that provides the highest PSNR and/or lowest nRMSD may be selected. In some such embodiments, if each of the PSNR and nRMSD provide different selections, and/or if either the PSNR or nRMSD
result in a tie, other factors may be considered in selecting a candidate. For example, one or more of the selection methods described below may be performed among some or all of the candidates, a candidate may be selected randomly, and/or a candidate that is "easier"
to use than other candidates (e.g. a candidate that is less computationally resource-intensive) may be selected. Illustrative examples of this approach are provided subsequently in the Example Application and Experimental Results sections.
[0110] In some embodiments, a signal profile is constructed based on specific signal characteristics so that comparisons can be made between the measured signal and the signal profile. By constructing the signal profile based on signal characteristics of particular interest, the comparisons may, in at least some embodiments, be biased towards those characteristics. For example, in the context of SCG, compressed sensing system 1 may provide and/or construct a signal profile which represents the morphology of a typical SCG signal (e.g. based on samples from one or more patients). A candidate may then be selected according to its similarity to the signal profile, as described below. Candidates with higher similarity to the signal profile tend to be more effective, in some circumstances, than candidates with lower similarity to the signal profile.
[0111] A signal profile reflects, encapsulates, or is otherwise based on known characteristics or morphologies of the particular signals being processed. For example, a signal profile for an SCG signal is likely to differ from a signal profile for an ECG signal, due to the different morphologies of those signals. A profile may be constructed statistically, e.g. by performing statistical time averaging. In such a construction, a sufficient number of samples of a particular signal of interest is acquired (this is sometimes referred to as the "training stage"). Then, windows or segments of the signals are time-averaged together, with appropriate timing synchronization or alignment (e.g. to align major peaks in the signal's amplitude). This averaging produces an overall template with a given window size. Such a construction may, for example, have two parameters: (1) the window size; and (2) the number of such windows used for averaging.
Preferably, the window size should be long enough to adequately characterize the signal of interest (e.g., in SCG one suitable window length may span at least 3 heart beats), but should not be too long so as to be affected by non-stationary distortions (e.g., in SCG, when under significant motion such as running, the signal's morphology may change significantly over 10-15 heart beats, and so it can be desirable to make the window size shorter than, for example, 10 heart beats). The number of windows may be based on statistical considerations. For example, the number of windows may be sufficiently large for the law of large numbers to apply, so that the averaging of the windows may bear a statistical relationship with the underlying distribution of the signal being measured.
For example, with SCG signals, an average of several hundred windows may provide a usable approximation.
[0112] Profiles may be constructed for an entire signal class (e.g., ECG class with a general PQRST profile) and/or a specific individual (e.g., an SCG of a particular individual with characteristic HIJK profile, for whom the ratios of the different wave amplitudes may be known). In some embodiments, multiple signal profiles may be generated for different conditions such as running, walking, sitting, etc. Of course, more specific characterizations may have more complex training requirements and/or computational or storage complexity, whereas less specific characterizations may be less accurate in some circumstances. For example, a physician may be interested in tracking the morphology of the R-wave for the ECG, or the H-wave for SCG. In such a case, a suitable profile could be created which includes windows with those specific characteristics. The windows may be collected over a sufficient number of patients so as to provide sufficient diversity (e.g., by processing signals from some or all of a physician's past or current patients).
[0113] In some embodiments where a signal profile is obtained, the signal profile may be used in the selection of a transformation matrix P. For example, a similarity metric based on the signal profile may be used instead of, or in addition to, the ranking matrices described above. Such similarity metrics include timing correlation metrics, which compare the timing of particular indicia in the signal profile with timing characteristics of the candidates. For example, in the particular case where the candidates are wavelets, a method with timing correlation has been described by B. N. Singh and A. K.
Tiwari, in the paper "Optimal selection of wavelet basis function applied to ECG signal denoising", Digit. Signal Process., vol. 16, no. 3, pp. 275-287, May 2006, which is hereby incorporated herein by reference for all purposes. In that paper, a timing correlation comparison is performed for the purpose of denoising. For example, revisiting the previous scenario of a physician interested in tracking the H-wave of an SCG
signal, a wavelet basis can be specifically selected for a specific user, based on the corresponding correlation values of the wavelet basis with the user's signal profile.
[0114] In some embodiments, a suitable basis (or bases) is constructed, rather than (or in addition to) being selected from a set of candidates. Such bases may, for example, be constructed using empirical decomposition and/or variants thereof, such as empirical mode decomposition, variational mode decomposition, and empirical wavelet transform (EWT).
EWT is described by J. Gilles in the paper "Empirical Wavelet Transform", IEEE
Trans.
Signal Process., vol. 61, no. 16, pp. 3999-4010, Aug. 2013, which is hereby incorporated herein by reference for all purposes. That paper aims to construct transformation matrices that incorporate knowledge of the modal nature of the measured signal. This approach evaluates the modal nature of the signal by measuring the number of local maxima in the signal. In the context of the present disclosure, the application of EWT may, in some embodiments, address two key parameters: (1) the placement of the modes; (2) the number of modes to be used. In some embodiments, the modes may be placed to target the spectrum morphology of an SCG signal. For example, in SCG, the modes may span the general SCG spectrum from 1Hz-100Hz. Modes above and below this range may be suppressed, as they may correspond to artifacts introduced by motion, baseline drift, equipment noise, and/or other extraneous interferences. In some embodiments, the number of modes used in the EWT construction is based on the desired number of basis functions.
For example, the mode number may be selected to accommodate the rank of the transformation matrix of a desired dimension. The above approach may allow for the construction of a transformation matrix that is specifically adapted to the defined signal profile, which may provide increased accuracy in reconstruction of the measured signal relative to matrix-selection approaches (such as those described above).
[0115] In block 165, compressed sensing system 1 obtains solution matrix 0. In some embodiments, 0 is obtained by finding the product 8 = . In such embodiments, decoder 50 is provided with at least one of sensing matrix (1) (with which solution matrix solution matrix e may be generated) and/or solution matrix e. As with sensing matrix (/) and transformation matrix tIJ, solution matrix 0 may be generated by decoder 50 after receiving encoded data y, as shown in Figure 2. In such cases, decoder 50 may use sensing matrix (/), which may be provided to it by encoder 30 or otherwise obtained (e.g. as part of initialization method 205). Alternatively, or in addition, solution matrix () may be generated prior to decoding method 150, e.g. in initialization method 205, as shown in Figure 3. In such cases, decoder 50 may not require access to sensing matrix 0. Solution matrix 0 may, for example, be generated before or during decoding, and may be generated by compressed sensing system 1 or by another system and/or user. A solution matrix 0 is provided or obtained for each (4),T) pair used to encode and/or decode data.
In some embodiments, compressed sensing system 1 is pre-programmed with one or more solution matrices 0 for use in decoding method 150.
[0116] In some embodiments where encoder 30 switches between different sensing matrices 0, encoder 30 may notify decoder 50 when a different sensing matrix will be used during encoding. Such notification may be provided, for example, along with encoded data y (e.g. as part of header data communicated with encoded data y) and/or separately (e.g. through a control channel).
[0117] Once the encoded signal y has been received and transformation matrix W
and solution matrix 0 have been obtained, decoder 50 goes to block 170 and finds a sparse solution to Eq. (6). That is, decoder 50 attempts to find the sparsest . that satisfies y = 0 , or an approximation thereof Recall that the optimal solution is guaranteed to be equal to s, as defined in Eq. (4), if 0 satisfies the restricted isometry property. If transformation matrix lIT and sensing matrix 0 are incoherent, as is likely in at least some of the above-disclosed embodiments, then 8 likely satisfies the restricted isometry property.
[0118] As noted above, finding the optimal (sparsest) solution is NP-hard, but reasonable approximations of the optimal solution can be found via numerical analysis. As will be understood by persons skilled in the art, several methods for solving such problems are available. These include, without limitation, basis pursuit, matching pursuit, and orthogonal matching pursuit. Techniques for performing at least some of these methods are disclosed by, for example, Y. Tsaig and D. L. Donoho in their 22 October 2004 paper entitled "Extensions of Compressed Sensing", and by J. A. Tropp and A. C.
Gilbert in their 11 April 2005 paper entitled "Signal Recovery From Partial Information Via Orthogonal Matching Pursuit"; both papers are incorporated herein by reference.
[0119] After finding a sparse solution , decoder 50 goes to block 175 to reconstruct the input signal. Since x = 11"s by construction, as shown in Eq. (4), decoder 50 may generate its output according to the equation:
i = Wg, (7) where is the sparse solution found in block 170 and W is the transformation matrix obtained in block 160.
[0120] In some embodiments, decoder 50 may receive from encoder 30 multiple encodings of a set of measurements. For example, encoder 30 may encode a set of measurements with a first sensing matrix 01 and may encode a copy of that set of measurements with a second sensing matrix 02; alternatively, or in addition, the different encodings may be of non-identical sets of measurements (e.g. the second set of encoded measurements may overlap with, or be a superset of subset of, the first set of encoded measurements). These multiple encodings of a set of measurements may be received by decoder 50, which may process some or all of them as indicated above at blocks 155-175.
Once each encoded set of measurements has been reconstructed to yield an approximation of the input signal, decoder 50 may optionally attempt to recover the input signal (and/or an approximation thereof) based on the multiple approximations at block 177.
For example, decoder 50 may perform an ensemble machine learning or committee machine approach, e.g., as described by S. K. Ambat, S. Chatterjee, and K. V. S. Hari, in their paper entitled "A Committee Machine Approach for Compressed Sensing Signal Reconstruction", IEEE Trans. Signal Process., vol. 62, no. 7, pp. 1705-1717, Apr. 2014, which is incorporated herein by reference for all purposes. Such approaches are likely to cause method 100 to require additional time, computational resources, and/or storage resources, but may improve the overall accuracy of method 100 in some circumstances. In particular, an ensemble learning approach may allow for a reduction in the amount of data that needs to be acquired or measured, but may (in some embodiments) require a more intricate receiver and/or decoder. For example, a compressed sensing system 1 which performs such an approach may measure fewer samples at the sensor 10 (thereby acquiring less data) and encode those measured samples in multiple ways at encoder 30.
These multiple encodings may be received at decoder 50, individually decoded, and then used at block 177 to recover an approximation which, in some circumstances, may be more accurate than one or more of the individually-decoded approximations produced at block 175.
[0121] After a solution has been reconstructed in block 175 (and, optionally, recovered from multiple reconstructions in block 177), decoder 50 then outputs the reconstructed signal 2 in block 180.2 may be further processed.
Example Application to SCG
[0122] As mentioned above, selection of a most appropriate transformation matrix W is strongly influenced by the morphology of the signal being processed. It is generally more difficult to find appropriate transformation matrices P for more complex signals (i.e.
signals that are less regular in their morphology), since it is more challenging to predict whether a given transformation matrix W will induce sparsity on such signals.
It is often expedient to conduct empirical research to determine whether certain transformation matrixes IP are likely to provide significant sparsity (and thus strong compression) for a given class of signal. Alternatively, or in addition, if prior knowledge about the signal morphologies are available, and if the required resources in storage and computation are available, then basis selection and/or construction methods may be performed, as described above.
[0123] One class of signals that presents particular challenges is SCG
signals. SCG is a physiological signal, similar in some respects to the electrocardiogram (ECG), the photoplethysmogram (PPG), or the electromyogram (EMG). Like the more conventional ECG, SCG measures signals produced by the human (or other animals') bodies that are related to heart activity. In particular, SCG measures the acceleration generated by the rhythmic contraction and relaxation of the heart on the chest wall.
[0124] Figure 4 shows SCG apparatus 40 according to an example embodiment.
Apparatus 40 includes a sensor 42 comprising an accelerometer 44.
Accelerometer 44 may be a three-axis accelerometer. Output signals from accelerometer 44 are conditioned by analog signal conditioning circuits 45X, 45Y and 45Z. Signal conditioning may comprise filtering and amplification for example. Conditioned signals are digitized by analog-to-digital converters (ADCs) 47X, 47Y and 47Z. In some embodiments a single ADC
may digitize two or more signals by sampling the signals at interleaved times.
[0125] Digitized signals are assembled into blocks of N samples in buffers 48X, 48Y and 48Z. Each block of N samples is then encoded using a suitable sensing matrix in blocks 49X, 49Y and 49Z which respectively apply sensing transformations (px, (py and (pz to yield reduced signals 50X, 50Y and 50Z (which may optionally be further compressed by data compression 52, e.g. via link-layer compression) and transmitted to a processor 55 by way of data link 54. Processor may apply data decompression 56 corresponding to data compression 52 to recover the reduced signals 50X, 50Y and 50Z. Processor 55 then applies a decoder 55A to process reduced signals 50X, 50Y and 50Z to recover the original signals digitized by ADCs 47 (or close approximations thereof).
Decoder 55A
finds a sparse solution at block 57 and reconstructs the signal at block 58 substantially as described above. These signals may be further processed, at block 59, according to any suitable SCG techniques (currently known or developed in the future) to yield the reconstructed signal.
[0126] An SCG accelerometer (i.e. the sensor) is usually placed on the center of the anterior chest wall, where the ribs come together. Signals measured by the sensor may be analyzed to determine, for example, the periodicity, strength, and coordination of contractions of the heart. SCG has attracted some attention for its usefulness in accurately identifying early coronary artery diseases, which are less readily identifiable in conventional ECGs or BCGs.
[0127] SCG records the compressional waves produced primarily by both heart wall motion and blood flow, resulting in a more data-rich signal than the conventional ECG. As a consequence, the ECG is generally much more sparse in the time domain than SCG, and transformations that may be applicable to the ECG or other physiological signals may have reduced efficacy in the context of SCG. In addition to the morphological complexity of SCG signals, it is not uncommon for sampling rates in SCG systems to be on the order of 1000 Hz or higher, as even very small and/or brief fluctuations in an SCG
signal can be useful in later analysis.
[0128] Although it is possible to perform SCG on a single axis (mono-axial SCG), it is sometimes preferable to perform SCG in three axes (tri-axial SCG), thereby providing a more complete view of the heart's activity (as measured from the position of the accelerometer). In tri-axial SCG, three axes are typically identified by convention: the x-axis points to the subject's right hand, the y-axis toward the subject's feet, and the z axis points outward from the subject's chest. The SCG signals acquired using these axes are referred to as x-SCG, y-SCG, z-SCG, respectively. Output of a typical tri-axial SCG is shown in Figures 4a, 4b, and 4c (corresponding to x-SCG, y-SCG, and z-SCG, respectively). A corresponding ECG is also plotted in Figure 4d, for reference. Although the ECG displays quasi-periodicity, the SCG signals tend to have greater variation. This is due to SCG accelerometers being sensitive to several sources of body mechanical movements, whereas the ECG depicts the relatively rhythmic electrical activities of the heart cycle. The example ECG signal in Figure 4d is highly regular, whereas the example signals from a tri-axial SCG sensor in Figures 4a to 4c are relatively quite irregular.
[0129] The large amounts of relatively sensitive data that may need to be processed in SCG systems may pose challenges to adoption of SCG in many circumstances.
These challenges are exacerbated, for example, in long-term care, remote sensing systems, and/or other situations where data is collected over long periods of time (increasing the amount of data to store/transmit/analyze) and/or transmitted over data communication channels where bandwidth concerns may be significant. For example, these challenges are often present when dealing with the long-term storage of SCG data, especially when that data corresponds to multiple patients. By compressing SCG data associated with patient records prior to storage and decompressing the data upon or after retrieval from storage, the burden on the data storage system can be reduced and an acceptable level of accuracy in the reconstructed signals can be maintained. As another example, in a telemedicine setting, transmission of compressed SCG data can help to conserve bandwidth and/or power, each of which may be scarce.
[0130] Although the discrete cosine transform (DCT) family of transformations is often used to construct basis functions in signal processing, conventionally-defined DCTs tend to produce large errors when used to decode SCG signals, even with modest compression ratios. In some embodiments, decoder 50 uses a modified "pseudo-DCT" (PDCT) to generate a transformation matrix IP. To construct a transformation matrix IP
using a PDCT, decoder 50 may construct an NxN matrix (where N is the size of the original vector x) with the kth column defined as:
kn-pdctk = N N
[1, cos (¨) cos ((N 1)klIT . (6) [0131] The columns of Eq. (6) are then normalized:
pdctk 1Pic = __ __ (7) and the mean of each of the columns (or, in some embodiments, each of the columns 2 to N) is set to zero if it is not already 0. In the PDCT as defined above, the mean of at least columns 2 to N is already 0 by construction, and so explicit removal of the columns' mean to 0 may not be necessary. Decoder 50 then sets W = ...
oPN]. This transformation matrix ir has a fixed basis and may thus be used with the i.i.d. Gaussian sensing matrix discussed above.
[0132] Other transformation matrices IP may be used to encode and decode SCG
signals.
For example, a transformation matrix W may be based on a wavelet transform (WT). In general, a WT is defined based on a generating function/wavelet CO. All other wavelet bases can be generated by translating and scaling the generating wavelet as follows:
1 (t¨b Oct,b(t) = , (8) VFL a where a is positive and defines the scale, and b defines the shift.
[0133] Not every WT basis is suitable for every type of signal. For SCG
signals in particular, we have determined that, among well-known WT bases, certain WTs produce reasonable results (and, in general, tend to better results than the conventional DCT), as quantified by the PSNR and nRMSD, defined below. These are the Haar, Daubechies 2, Daubechies 5, Symlet 2, Symlet 4, Biorthogonal 4.4, Biorthogonal 5.5, Reverse Biorthogonal 4.4, and Reverse-Biorthogonal 5.5 wavelet transforms.
Compression Ratios [0134] The above transforms are suitable for both mono-axial and tri-axial SCG. Due to the relatively complicated form of the SCG, even applications involving mono-axial SCG
generate large amounts of data and may be challenging in terms of storage, transmission and computational complexity. Tr-axial SCG is particularly data-intensive (as shown in Figure 4), and may therefore benefit greatly from the application of the above methods and systems.
[0135] In experimental testing, compression ratios as low as 'A (i.e. M ===.; -41 N, where M is the size of the encoded vector y and N is the size of the input vector x) have been reached without unacceptably degrading the accuracy of the reconstructed signal. Even more aggressive compression (e.g. a compression ratio of 1/5) has yielded acceptable results on some input vectors, and particularly on low-frequency input vectors. In at least some of the tests conducted, each stream was compressed and decompressed with the same sensing and transformation matrices.
[0136] In the case of tri-axial SCG, the three-dimensional sensor readings are often split into three data streams (i.e. three different sets of input vectors), corresponding to x-SCG, y-SCG, and z-SCG. Signals in some axes may have different frequencies than signals in other axes; for example, in some circumstances the y-SCG axis may tend to have a higher-frequency signal and/or may be less amenable to aggressive compression. In these or other situations, a different (e.g. less aggressive, larger) compression ratio may be used for the y-SCG stream than for one or both of the x-SCG and z-SCG streams. Similarly, because the z-SCG axis is often of particular interest to clinicians, a less aggressive (i.e. lower) compression ratio may be used for that stream to provide greater accuracy. The variation of compression ratios for one or more axes may vary depending on the circumstances. For example, the orientation of the patient (e.g. lying down vs. standing up) may influence the compression ratio selected for a particular axis. As a further example, certain diseases may result in indicia that are more prominent in certain axes, which may make increased accuracy in those axes more desirable (and therefore may affect the selection of compression ratios in those and/or other axes).
[0137] In the case of tri-axial SCG, one or more of the data streams may be encoded using constructed bases, e.g. as described above (for instance, by using EWT). For example, the y-SCG and z-SCG streams may use constructed bases, whereas the x-SCG stream may use a selected basis. Since it is often less resource-intensive (in terms of computational power) to select bases than to construct them, in some embodiments only high-priority streams (e.g. those of particular clinical significance) use constructed bases, whereas secondary or ancillary streams use selected bases for the sake of efficiency. Additionally, or alternatively, when different streams each use different constructed bases, those bases may have different numbers of modes. For example, in an application where the y-SCG stream is of principal clinical significance and the z-SCG stream is less clinically significant, the y-SCG stream may use a constructed basis with a greater number of modes (and thus, in some circumstances, may have greater accuracy in the encoding and decoding process) than the z-SCG stream.
101381 Compression ratios can be set by generating an M x N sensing matrix 4).
The ratio M/N is equal to the compression ratio. N is determined according to the size of the input vector (e.g. which may be set according to the size of the windows used to segment the input stream), but M may be chosen arbitrarily, subject to the proviso that aggressive compression may adversely affect the accuracy of the reconstructed signal. In embodiments where sensing matrix 4) is generated randomly (or pseudo-randomly), generation of a new sensing matrix 4) can be done relatively easily. The same N x N
transformation matrix may be used with the new sensing matrix 0.
Alternatively, or in addition, a new N x N transformation matrix W may be used, but in either case a new solution matrix 0 must be obtained.
101391New sensing matrices 4) may be selected based on the desired compression ratio ¨
e.g. if compressed sensing system 1 and/or a user believe that significant compression is possible without unacceptably impacting the accuracy of the decoded signal, the new sensing matrix 4) may be generated with a small M. Since the lower bound for Mhos a lower bound of K log NiK, the selection of M will, in some cases, be at least partially based on a prediction of the sparseness of x under transformation 41. In other words, the better-matched that transformation matrix tY is to the input signal, the smaller that M can be before unacceptably affecting decoding accuracy. Thus, the size of M (and therefore the size of 4)) is indirectly connected to the morphology of x.Preprocessing (Windowing).
[0140] When implementing compressed sensing as described above, e.g. in a tri-axial SCG system, the "raw" SCG signals on each axis (e.g. x-SCG, y-SCG, and z-SCG) may be segmented according to a windowing scheme (e.g. processing 1000 samples collected over a one-second window). The above methods do not generally require significant preprocessing of signal readings prior to encoding. Generally, larger window sizes (resulting in larger input vectors) result in increased computational complexity and, CA 02851028 2014-05-02 .
according to experimentation, slightly increased accuracy in some cases.
[0141] In experimental testing, using a window length of at least one complete SCG cycle (i.e. a complete heartbeat) provided adequate compressibility. In general, window lengths may be selected to be sufficiently long to encapsulate the complete SCG
morphologies, but not so long that non-stationary statistical behavior prevents accurate signal estimation.
Longer window cycles (e.g two complete SCG cycles) may produce better compression results, but may also result in increased computational complexity without, in some circumstances, a commensurate improvement in overall accuracy.
Experimental Results [0142] There are several methods available for measuring the results of encoding schemes generally. The present disclosure has been tested using at least two metrics, namely peak signal-to-noise ratio (PSNR) and normalized root-mean-square deviation (nRMSD). PSNR
may be derived as follows:
MAX
PSNR(dB)= 20 log Ihrfs1.., (9) where, for a group of signals X and their recoveries i, MAX is the maximum value of X
and MSE is the mean square error, calculated as follows:

MSE = ¨NTII[X(i,j) ¨ -X(i,j)12 .(10) i=o j=o [0143] A higher PSNR generally indicates better reconstruction performance.
The nRMSD
may be calculated according to the following formula:
nRMSD(%) = 100i ¨1 N
X - i NI[max(x) ¨ min(x) , (11) where x is an input signal and 2 is its recovered signal and both are vectors with size N x 1. The nRMSD is usually expressed in percentage, with a lower value generally indicating better system performance.
[0144] In one experiment, sample signals were taken from a dataset of SCG
signals. The sample signals were encoded and decoded using the compressed sensing system described above, using each of the above-identified fixed bases for the transformation matrix W. A
range of compression ratios were used in order to test the relative quality of each fixed basis.
[0145] The results of the experiment where the compression ratio was set to 1 (i.e. N = M, resulting in no compression) are shown in Figure 5a. Figure 5a shows three graphs, corresponding to (from top to bottom) PSNR results for x-SCG, y-SCG, and z-SCG.
Figure 5b is similarly arranged and shows the results of a similar experiment where the compression ratio was set to 1/6 (i.e. M = 1/6 N), which is typically a very aggressive compression ratio. Figures 5a and 5b show a poor performance by the conventional DCT, and relatively more robust performances by the PDCT and the specific WT bases that are discussed above.
[0146] PSNR was measured for sample data with several compression ratios in addition to the above, providing the results shown in Figure 6a. nRMSD was used to measure results of compression of the same input data at a range of compression ratios. Since nRMSD
works on individual signals (as opposed to groups of signals), these results were reported as an average across all of the sample signals (Figure 6b), the maximum results (Figure 6c), and the minimum results (Figure 6d). The results for conventional DCT are omitted from Figures 6a-d, as they dramatically skewed the graphs.
[0147] Although the results varied between the fixed bases, the obtained results show that, in certain circumstances, the above-described compressed-sensing system can achieve reasonably accurate results at compression ratios of 1/4 or even lower, with some bases providing acceptable results with an aggressive compression ratio of 1/6.
Interpretation of Terms [0148] Unless the context clearly requires otherwise, throughout the description and the claims:
= "comprise", "comprising", and the like are to be construed in an inclusive sense, as opposed to an exclusive or exhaustive sense; that is to say, in the sense of "including, but not limited to";
= "connected", "coupled", or any variant thereof, means any connection or coupling, either direct or indirect, between two or more elements; the coupling or connection between the elements can be physical, logical, or a combination thereof;
= "herein", "above", "below", and words of similar import, when used to describe this specification, shall refer to this specification as a whole, and not to any particular portions of this specification;
= "or", in reference to a list of two or more items, covers all of the following interpretations of the word: any of the items in the list, all of the items in the list, and any combination of the items in the list;
= the singular forms "a", "an", and "the" also include the meaning of any appropriate plural forms.
[0149] Words that indicate directions such as "vertical", "transverse", "horizontal", "upward", "downward", "forward", "backward", "inward", "outward", "vertical", "transverse", "left", "right", "front", "back", "top", "bottom", "below", "above", "under", and the like, used in this description and any accompanying claims (where present), depend on the specific orientation of the apparatus described and illustrated.
The subject matter described herein may assume various alternative orientations.
Accordingly, these directional terms are not strictly defined and should not be interpreted narrowly.
[0150] Although the above description refers at various times to vectors, persons skilled in the art will appreciate that the above methods and systems can be adapted to operate on matrices as well. For example, the above methods can be applied column-wise to a matrix of input sensor readings.

[0151] Embodiments of the invention may be implemented using specifically designed hardware, configurable hardware, programmable data processors configured by the provision of software (which may optionally comprise "firmware") capable of executing on the data processors, special purpose computers or data processors that are specifically programmed, configured, or constructed to perform one or more steps in a method as explained in detail herein and/or combinations of two or more of these.
Examples of specifically designed hardware are: logic circuits, application-specific integrated circuits ("ASICs"), large scale integrated circuits ("LSIs"), very large scale integrated circuits ("VLSIs"), and the like. Examples of configurable hardware are: one or more programmable logic devices such as programmable array logic ("PALs"), programmable logic arrays ("PLAs"), and field programmable gate arrays ("FPGAs")). Examples of programmable data processors are: microprocessors, digital signal processors ("DSPs"), embedded processors, graphics processors, math co-processors, general purpose computers, server computers, cloud computers, mainframe computers, computer workstations, and the like. For example, one or more data processors in a control circuit for a device may implement methods as described herein by executing software instructions in a program memory accessible to the processors.
[0152] Processing may be centralized or distributed. Where processing is distributed, information including software and/or data may be kept centrally or distributed. Such information may be exchanged between different functional units by way of a communications network, such as a Local Area Network (LAN), Wide Area Network (WAN), or the Internet, wired or wireless data links, electromagnetic signals, or other data communication channel.
[0153] For example, while processes or blocks are presented in a given order, alternative examples may perform routines having steps, or employ systems having blocks, in a different order, and some processes or blocks may be deleted, moved, added, subdivided, combined, and/or modified to provide alternative or subcombinations. Each of these processes or blocks may be implemented in a variety of different ways. Also, while processes or blocks are at times shown as being performed in series, these processes or blocks may instead be performed in parallel, or may be performed at different times.

[0154] The invention may also be provided in the form of a program product.
The program product may comprise any non-transitory medium which carries a set of computer-readable instructions which, when executed by a data processor, cause the data processor to execute a method of the invention. Program products according to the invention may be in any of a wide variety of forms. The program product may comprise, for example, non-transitory media such as magnetic data storage media including floppy diskettes, hard disk drives, optical data storage media including CD ROMs, DVDs, electronic data storage media including ROMs, flash RAM, EPROMs, hardwired or preprogrammed chips (e.g., EEPROM semiconductor chips), nanotechnology memory, or the like. The computer-readable signals on the program product may optionally be compressed or encrypted.
[0155] In some embodiments, the invention may be implemented in software. For greater clarity, "software" includes any instructions executed on a processor, and may include (but is not limited to) firmware, resident software, microcode, and the like. Both processing hardware and software may be centralized or distributed (or a combination thereof), in whole or in part, as known to those skilled in the art. For example, software and other modules may be accessible via local memory, via a network, via a browser or other application in a distributed computing context, or via other means suitable for the purposes described above.
[0156] Where a component (e.g. a software module, processor, assembly, device, circuit, etc.) is referred to above, unless otherwise indicated, reference to that component (including a reference to a "means") should be interpreted as including as equivalents of that component any component which performs the function of the described component (i.e., that is functionally equivalent), including components which are not structurally equivalent to the disclosed structure which performs the function in the illustrated exemplary embodiments of the invention.
[0157] Specific examples of systems, methods and apparatus have been described herein for purposes of illustration. These are only examples. The technology provided herein can be applied to systems other than the example systems described above. Many alterations, modifications, additions, omissions, and permutations are possible within the practice of this invention. This invention includes variations on described embodiments that would be apparent to the skilled addressee, including variations obtained by: replacing features, elements and/or acts with equivalent features, elements and/or acts; mixing and matching of features, elements and/or acts from different embodiments; combining features, elements and/or acts from embodiments as described herein with features, elements and/or acts of other technology; and/or omitting combining features, elements and/or acts from described embodiments.
[0158] It is therefore intended that the following appended claims and claims hereafter introduced are interpreted to include all such modifications, permutations, additions, omissions, and sub-combinations as may reasonably be inferred. The scope of the claims should not be limited by the preferred embodiments set forth in the examples, but should be given the broadest interpretation consistent with the description as a whole.

Claims (62)

WHAT IS CLAIMED IS:

1. A method of performing compressed sensing seismocardiography, the method comprising:
receiving an input set of data representing a seismocardiogram signal;
obtaining a sensing matrix;
determining an encoded set of data based on a product of at least a portion of the input set of data and the sensing matrix;
obtaining a transformation matrix corresponding to the sensing matrix;
obtaining a solution matrix based on a product of the sensing matrix and the transformation matrix, the solution matrix satisfying the restricted isometry property;
determining a sparse interim set of data so that a product of the sparse interim set of data and the solution matrix corresponds to the encoded set of data;
and determining an output set of data based on a product of the transformation matrix and the sparse solution so that the output set of data approximates the input set of data.

2. A method according to claim 1 wherein the sensing matrix has a size M
× N where N is a size of the input set of data and M is less than N.

3. A method according to any one of claims 1 and 2 wherein the sensing matrix is constructed from independent and identically distributed Gaussian distributions.

4. A method according to any one of claims 1 to 3 wherein obtaining the sensing matrix comprises retrieving the sensing matrix from a memory.

5. A method according to any one of claims 1 to 3 wherein obtaining the sensing matrix comprises generating the sensing matrix.

6. A method according to any one of claims 1 to 3 wherein obtaining the sensing matrix comprises receiving the sensing matrix from a user.

7. A method according to any one of claims 1 to 6 wherein the transformation matrix has a size N × N where N is a size of the input set of data.

8. A method according to any one of claims 1 to 7 wherein the transformation matrix and the sensing matrix are incoherent.

9. A method according to any one of claims 1 to 8 wherein the transformation matrix has a fixed basis.

10. A method according to claim 9 wherein the fixed basis comprises a PDCT
so that the k th column of the fixed basis is based on the construction given by:
where N is a size of the input set of data.

11. A method according to claim 10 wherein the fixed basis has normalized zero-mean columns.

12. A method according to claim 9 wherein the fixed basis comprises the Haar wavelet transform.

13. A method according to claim 9 wherein the fixed basis comprises the Daubechies 2 wavelet transform.

14. A method according to claim 9 wherein the fixed basis comprises the Daubechies 5 wavelet transform.

15. A method according to claim 9 wherein the fixed basis comprises the Symlet 2 wavelet transform.

16. A method according to claim 9 wherein the fixed basis comprises the Symlet 4 wavelet transform.

17. A method according to claim 9 wherein the fixed basis comprises the Bioorthogonal 4.4 wavelet transform.

18. A method according to claim 9 wherein the fixed basis comprises the Bioorthogonal 5.5 wavelet transform.

19. A method according to claim 9 wherein the fixed basis comprises the Reverse Bioorthogonal 4.4 wavelet transform.

20. A method according to claim 9 wherein the fixed basis comprises the Reverse Bioorthogonal 5.5 wavelet transform.

21. A method according to any one of claims 1 to 20 wherein obtaining the transformation matrix comprises retrieving the transformation matrix from a memory.

22. A method according to any one of claims 1 to 20 wherein obtaining the transformation matrix comprises generating the transformation matrix.

23. A method according to any one of claims 1 to 20 wherein obtaining the transformation matrix comprises receiving the transformation matrix from a user.

24. A method according to any one of claims 1 to 23 wherein determining a sparse interim set of data comprises finding a sparse solution to:
where y is the encoded set of data and .PHI. is the sensing matrix.

25. A method according to claim 21 wherein finding a sparse solution comprises performing basis pursuit.

26. A method according to claim 21 wherein finding a sparse solution comprises performing matching pursuit.

27. A method according to claim 21 wherein finding a sparse solution comprises performing orthogonal matching pursuit.

28. A method according to any one of claims 1 to 27 wherein the seismocardiogram signal is a tri-axial seismocardiogram signal.

29. A method according to any one of claims 1 to 28 wherein the input set of data is segmented according to a windowing scheme.

30. A method according to any one of claims 1 to 29 comprising:
obtaining a second sensing matrix;
using the second sensing matrix in place of the sensing matrix;
obtaining a second solution matrix based on the second sensing matrix, the second solution matrix satisfying the restricted isometry property; and using the second solution matrix in place of the solution matrix.

31. A method according to claim 30 comprising determining an accuracy measure of the output set of data, wherein the second sensing matrix is obtained in response to the determination of the accuracy measure.

32. A method according to any one of claims 30 and 31 wherein the second sensing matrix is obtained in response to user input.

33. A method according to any one of claims 30 to 32 wherein the second sensing matrix is obtained based on a desired compression ratio.

34. A system for compressed sensing seismocardiography, the system comprising:
an encoder for encoding an input set of data, the encoder configured to determine an encoded set of data based on a product of at least a portion of the input set of data and a sensing matrix; and a decoder for decoding the encoded set of data, the decoder configured to:
determine a sparse interim set of data so that a product of the sparse interim set of data and a solution matrix corresponds to the encoded set of data; and determine an output set of data based on a product of a transformation matrix and the sparse solution so that the output set of data approximates the input set of data;

wherein:
the transformation matrix corresponds to the sensing matrix;
the solution matrix is based on a product of the sensing matrix and the transformation matrix; and the solution matrix satisfies the restricted isometry property.

35. A system according to claim 34 wherein the sensing matrix has a size M
× N
where N is a size of the input set of data and M is less than N.

36. A system according to any one of claims 31 and 35 wherein the sensing matrix is constructed from independent and identically distributed Gaussian distributions.

37. A system according to any one of claims 31 to 36 wherein the transformation matrix has a size N × N where N is a size of the input set of data.

38. A system according to any one of claims 31 to 37 wherein the transformation matrix and the sensing matrix are incoherent.

39. A system according to any one of claims 31 to 38 wherein the transformation matrix has a fixed basis.

40. A system according to claim 39 wherein the fixed basis comprises a PDCT
so that the k th column of the fixed basis is based on the construction given by:
where N is a size of the input set of data.

41. A system according to claim 40 wherein the fixed basis has normalized zero-mean columns.

42. A system according to claim 39 wherein the fixed basis comprises the Haar wavelet transform.

43. A system according to claim 39 wherein the fixed basis comprises the Daubechies 2 wavelet transform.

44. A system according to claim 39 wherein the fixed basis comprises the Daubechies wavelet transform.

45. A system according to claim 39 wherein the fixed basis comprises the Symlet 2 wavelet transform.

46. A system according to claim 39 wherein the fixed basis comprises the Symlet 4 wavelet transform.

47. A system according to claim 39 wherein the fixed basis comprises the Bioorthogonal 4.4 wavelet transform.

48. A system according to claim 39 wherein the fixed basis comprises the Bioorthogonal 5.5 wavelet transform.

49. A system according to claim 39 wherein the fixed basis comprises the Reverse Bioorthogonal 4.4 wavelet transform.

50. A system according to claim 39 wherein the fixed basis comprises the Reverse Bioorthogonal 5.5 wavelet transform.

51. A system according to any one of claims 31 to 50 wherein to determine a sparse interim set of data comprises finding a sparse solution to:
where y is the encoded set of data and .PHI. is the sensing matrix.

52. A system according to claim 51 wherein finding a sparse solution comprises performing basis pursuit.

53. A system according to claim 51 wherein finding a sparse solution comprises performing matching pursuit.

54. A system according to claim 51 wherein finding a sparse solution comprises performing orthogonal matching pursuit.

55. A system according to any one of claims 31 to 54 wherein the seismocardiogram signal is a tri-axial seismocardiogram signal.

56. A system according to any one of claims 31 to 55 wherein the input set of data is segmented according to a windowing scheme.

57. A system according to any one of claims 31 to 56, wherein the system is configured to:
obtain a second sensing matrix;
use the second sensing matrix in place of the sensing matrix;
obtain a second solution matrix based on the second sensing matrix, the second solution matrix satisfying the restricted isometry property; and use the second solution matrix in place of the solution matrix.

58. A system according to claim 57 wherein the system is configured to determine an accuracy measure of the output set of data and to obtain the second sensing matrix in response to the determination of the accuracy measure.

59. A system according to any one of claims 57 and 58 wherein the system is configured to obtain the second sensing matrix in response to user input.

60. A system according to any one of claims 57 to 59 wherein the system is configured to obtain the second sensing matrix based on a desired compression ratio.

61. An apparatus having any new and inventive feature, combination of features, or sub-combination of features as described herein.

62. A method having any new and inventive steps, acts, combination of steps and/or acts or sub-combination of steps and/or acts as described herein.