Related Applications

[0001] This application claims the benefit of the priority of US application No. 61/837547 filed 20 June 2013, US application No. 61/954533 filed 17 March 2014 and US application No. 14/247070 filed 7 April 2014, all of which are incorporated herein by reference.

Technical Field

[0002] The subject matter of this disclosure relates to estimating a change in an inhalation- exhalation ratio (IER) of human or other animal subjects. Particular aspects of the invention provide systems and methods for using sensing systems to sense and/or monitor motion of such subjects' lungs from which the change in IER may be determined.

Background

[0003] There is a general desire to use sensing systems (such as Impulse-Radio Ultra Wideband (IR-UWB or, for brevity, UWB) radar systems or Ballistocardiography (BCG) systems) to sense and/or monitor the vital signs of human and/or other animal subjects to detect changes in IER. One non-limiting reason for this desire is that existing techniques (such as auscultation and spirometry with a pneumotachograph) may require contact with and/or participation of the subject, which can be undesirable in some clinical contexts. As another example, it can be inconvenient to conduct existing techniques over long periods (e.g. days or weeks) due to active involvement of the clinician and/or subject. Existing techniques may be inaccurate, costly, inconvenient to apply, and/or may involve devices which require frequent calibration.

[0004] A hospital doctor treating a recently admitted asthma attack victim typically wants to know if an ordered therapy is lessening the severity of asthma symptoms. Often these patients are not in the intensive care unit, and therefore are not receiving continuous monitoring. Instead, a medical professional visits the patient periodically (typically once every one or more hours) and makes an assessment based on a number of objective measurements (none of which (alone) indicate a change in symptoms), such as heart and respiration rates and oxygen saturation, and a subjective observation of the work that the patient is performing to breathe, which is based on use of accessory muscles and chest wall recession. These infrequent and episodic measurements may miss trends because of natural short term variations, making therapy assessment problematic. At home, a caregiver of an asthmatic child often first becomes aware of an impending attack by the sound of wheezing. However by the time the wheezing is audible, it is usually too late to avert the attack and may be too late to avoid an emergency department visit. In addition, an asthma attack can progress without audible wheezing, and an audible wheeze may not occur until asthma treatment and improved airflow occurs.

[0005] An indicator that is well known to con-elate well with worsening (improving) asthma symptoms is a decrease (increase) in the IER. Both of the foregoing examples are exemplary of the general desire to use sensing systems to detect changes in IER.

[0006] The foregoing examples of the related art and limitations related thereto are intended to be illustrative and not exclusive. Other limitations of the related art will become apparent to those of skill in the art upon a reading of the specification and a study of the drawings.

Summary

[0007] The following embodiments and aspects thereof are described and illustrated in conjunction with systems, tools and methods which are meant to be exemplary and illustrative, not limiting in scope. In various embodiments, one or more of the above- described problems have been reduced or eliminated, while other embodiments are directed to other improvements.

[0008] One aspect of the invention provides a method for estimating a change in an inhalation-exhalation ratio (IER) of a subject based on a reflected signal received by a UWB radar system. The method comprises processing a signal received from the sensing system to obtain a plurality of frequency domain representations of the signal. The method also comprises discerning, for each frequency domain representation, a fundamental characteristic of the frequency domain representation corresponding to a fundamental frequency corresponding to motion of one or more of the subject's lungs and at least one harmonic characteristic of the frequency domain representation corresponding to at least one harmonic frequency corresponding to the motion of one or more of the subject's lungs. The method further comprises estimating the change in the IER based on one or more fundamental characteristics and one or more harmonic characteristics of the plurality of frequency domain representations.

[0009] In some embodiments of this or other aspects, processing the signal received from the sensing system comprises determining an estimated power spectral density of the signal and at least one of the plurality of frequency domain representations of the signal comprises the determined power spectral density estimate.

[0010] In some embodiments of this or other aspects, estimating the change in the IER based on one or more fundamental characteristics and one or more harmonic characteristics comprises estimating the change in the IER based on a change in one or more ratios, each of the one or more ratios based on the fundamental characteristic and the at least one harmonic characteristics of at least one of the plurality of frequency domain representations. In some embodiments, at least one of the one or more ratios comprises a ratio average. The ratio average may be based on a first ratio based on a fundamental characteristic and a harmonic characteristic of a first frequency domain representation and a second ratio based on a fundamental characteristic and a harmonic characteristic of a second frequency domain representation. In some embodiments, the ratio average comprises a moving average of ratios of fundamental characteristics and harmonic characteristics of the plurality of frequency domain representations. In some embodiments, the ratio average comprises a linear regression of ratios of fundamental characteristics and harmonic characteristics of the plurality of frequency domain representations.

[0011] In some embodiments of this or other aspects, the plurality of frequency domain representations are based on a plurality of temporal windows of the signal. Each of the plurality of temporal windows may be associated with a time domain representation of the signal and at least one of the plurality of frequency domain representations of the signal. A . first time domain representation of the signal associated with a first temporal window of the plurality of temporal windows may overlap with a second time domain representation of the signal associated with a second temporal window of the plurality of temporal windows. At least one of the plurality of temporal windows may comprise a Hamming window.

[0012] In some embodiments of these or other aspects, processing the signal may comprise determining a moving average based on one or more windows of the plurality of windows. Obtaining a frequency domain representation of the signal may comprise omitting data from the frequency domain representation corresponding to a portion of the signal having a magnitude that differs from the moving average by more than a moving threshold value.

[0013] In some embodiments of this or other aspects, discerning, for each of the frequency domain representations of the signal, the fundamental characteristic and the at least one harmonic characteristics comprises identifying one or more peaks in the frequency domain representation.

[0014] In some embodiments of this or other aspects, the sensing system comprises a UWB radar system for directing UWB pulses toward a subject, receiving reflected pulses, and generating a signal based on the reflected pulses. In some embodiments of this or other aspects, the sensing system comprises a BCG sensor.

[0015] In some embodiments of this or other aspects, estimating the change in the IER based on one or more fundamental characteristics and one or more harmonic characteristics comprises determining a first change based on one or more fundamental characteristics and one or more harmonic characteristics. Estimating the change in the IER may also comprise determining a second change based on one or more fundamental characteristics and one or more harmonic characteristics. Estimating a change in the IER may further comprise determining that the first change is less than a change threshold and that the second change is greater than the change threshold. Estimating a change in the IER may still further comprise estimating the change in IER based on the first change and omitting the second change from the estimation of the change in IER.

[0016] Some embodiments may determine that the first change is less than the change threshold and greater than a warning threshold. Such embodiments may provide to a user a warning based on the first change. Some embodiments may determine a substitute change based on the second change. The substitute change may omit at least one of the one or more fundamental characteristics and the one or more harmonic characteristics on which the second change is based. Such embodiments may estimate the change in the IER based on the substitute change.

[0017] In some embodiments of this or other aspects, at least one of the fundamental characteristic and the at least one harmonic characteristic of a frequency domain

representation comprises a power of a peak of the frequency domain representation. The peak may correspond to at least one of the fundamental frequency and the at least one harmonic frequencies. The power of the peak may be determined based on a sum of values of such frequency domain representation in a vicinity of at least one of the fundamental frequency and the at least one harmonic frequencies. [0018] In some embodiments of this or other aspects, at least one of the fundamental characteristic and the at least one harmonic characteristic of a frequency domain

representation comprises a power of the frequency domain representation in a neighbourhood of the frequency domain representation. The neighbourhood may correspond to at least one of the fundamental frequency and the at least one harmonic frequencies. The power of the neighbourhood may be determined based on a sum of values of such frequency domain representation in the neighbourhood.

[0019] In some embodiments of this or other aspects, respiration rate information is received from an external system. Discerning, for each of the frequency domain representations of the signal, the fundamental characteristic and the at least one harmonic characteristic may be based at least in part on the received respiration rate information.

[0020] In some embodiments of this or other aspects, respiration rate information is received from an external system. For each of the frequency domain representations of the signal, the fundamental frequency and the at least one harmonic frequency may be estimated based at least in part on the received respiration rate information.

[0021] Another aspect of the invention provides a system for estimating a change in an inhalation-exhalation ratio (IER) of a subject. The system comprises a sensing system and a processor (which may be one or more processors) connected to receive a signal from the sensing system. The processor is configured to process a signal received from the sensing system to obtain a plurality of frequency domain representations of the signal. The processor is also configured to discern, for each frequency domain representation, a fundamental characteristic of the frequency domain representation corresponding to a fundamental frequency corresponding to motion of one or more of the subject's lungs and at least one harmonic characteristic of the frequency domain representation corresponding to at least one harmonic frequency corresponding to the motion of one or more of the subject's lungs. The processor is further configured to estimate the change in the IER based on one or more fundamental characteristics and one or more harmonic characteristics of the plurality of frequency domain representations.

[0022] Systems according to some embodiments may comprise one or more processors configured to perform any of the methods described herein. In some embodiments, the sensing system comprises a UWB radar system for directing UWB pulses toward a subject, receiving reflected pulses, and generating a signal based on the reflected pulses, in some embodiments, the sensing system comprises a BCG sensing system.

[0023] Another aspect of the invention provides a computer readable medium containing program instructions for estimating a change in an inhalation-exhalation ratio (IER) of a subject based on a reflected signal received by a sensing system, wherein execution of the program instructions by one or more processors of a computer system causes the one or more processors to carry out the steps of a method. The method comprises processing a signal received from the sensing system to obtain a plurality of frequency domain representations of the signal. The method also comprises discerning, for each frequency domain representation, a fundamental characteristic of the frequency domain representation corresponding to a fundamental frequency corresponding to motion of one or more of the subject's lungs and at least one harmonic characteristic of the frequency domain representation corresponding to at least one harmonic frequency corresponding to the motion of one or more of the subject's lungs. The method further comprises estimating the change in the IER based on one or more fundamental characteristics and one or more harmonic characteristics of the plurality of frequency domain representations.

[0024] In addition to the exemplary aspects and embodiments described above, further aspects and embodiments will become apparent by reference to the drawings and by study of the following detailed descriptions.

Brief Description of the Drawings

[0025] Exemplary embodiments are illustrated in referenced figures of the drawings. It is intended that the embodiments and figures disclosed herein are to be considered illustrative rather than restrictive.

[0026] Figure 1 A is a block diagram of an ultra-wideband (UWB) sensing system for estimating one or more physiological characteristics of a human or other animal according to a particular exemplary embodiment. Figure IB is a schematic view showing the Figure 1 A UWB sensing system estimating one or more physiological characteristics of a human subject through a mattress by way of UWB pulses.

[0027] Figure 2A is a block diagram representation of a method for estimating a change in IER of a human or other animal according to a particular embodiment which may be performed by the Figure 1A UWB sensing system. [0028] Figure 2B is a schematic depiction of a filter that may be applied to data received by the radar system of the Figure 1A sensing system for separating the respiration fundamental from the heart rate fundamental of a subject and for removing DC components.

[0029] Figure 3 shows an example of the type of frequency domain data that could be obtained as a part of the Figure 2A estimation method according to some embodiments.

[0030] Figure 4A shows an example of the variation of the magnitude of the fundamental with IER and the corresponding, opposing variation of the normalized magnitude of the harmonic characteristic with IER using example frequency domain data derived from example measurements of the Figure 1A UWB sensing system.

[0031] Figure 4B shows an example of the variation of the magnitude of the fundamental with IER and the corresponding, opposing variation of the normalized magnitude of the harmonic characteristic with IER using example frequency domain data derived from example measurements of an example BCG sensing system.

[0032] Figure 5A is a block diagram representation of a method for obtaining frequency domain data that may be used in the Figure 2A method for estimating a change in IER according to some embodiments.

[0033] Figure 5B is a schematic depiction of a technique for parsing time domain measurement data into temporal windows, transforming the data from the time domain windows into the frequency domain and averaging the windows of the corresponding frequency domain data which may be used in the method of Figure 5 A according to a particular embodiment.

[0034] Figure 6 is a table of example data set of ratios and their corresponding ratio averages which may be used in the Figure 2A method of estimating a change in IER according to an example embodiment.

Detailed Description

[0035] Throughout the following description specific details are set forth 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. Accordingly, the description and drawings are to be regarded in an illustrative, rather than a restrictive, sense.

[0036] One technique for early detection of lower airway obstruction (LAO) symptoms is measurement of the subject's inhalation-expiration ratio (IER). IER is the ratio of inspiratory respiration time (¾) to expiratory respiration time (¾), which typically ranges from 1 :2 to 1 :3 in healthy human children. Since prolonged exhalation (relative to inhalation) is a clinical sign commonly used to detect LAO symptoms, a low IER ratio outside of a healthy range for a given subject can be a useful quantitative indicator of LAO symptoms.

[0037] There is a general desire to use sensing systems to sense and/or monitor the vital signs of human and/or other animal subjects to detect changes in IER. One example of such a sensing system is Impulse-Radio Ultra Wideband (IR-UWB or, for brevity, UWB) radar. One non-limiting reason that IR-UWB is desirable is that it can be configured to not interfere with (or be interfered by) other wireless communication equipment. Also, or alternatively, IR- UWB radar can be used without contacting or otherwise inconveniencing the subject, does not require subject participation, and can be used at relatively low power, suitable for continuous monitoring. IR-UWB sensing systems may, for example, be placed under a subject's bed sheets, in the back or head rest of a diagnostic chair or similar device, on a ceiling above a subject's bed, and/or otherwise positioned on or near a subject for long-term and/or continuous monitoring. Another example of such a sensing system is the

ballistocardiogram (BCG), which may be placed in contact with a subject, for example by placement in the back or head rest of a diagnostic chair or similar device, and/or by mechanical coupling to the patient. The BCG may enable continuous and/or longer-term measurement without necessarily inconveniencing the subject or requiring subject participation.

[0038] Aspects of the invention provide systems and methods for estimating a change in an inhalation-exhalation ratio (IER) of a subject based on a signal received by a sensing system. Time domain measurement data may be obtained using a UWB radar system, BCG sensing system, and/or another sensing system. The time domain measurement data may be converted to the frequency domain to obtain a plurality of frequency domain representations of the measurement data (or, for brevity, frequency domain data). The system may comprise a processor configured to (or the method may comprise) process the frequency domain data to discern characteristics of the frequency domain data corresponding to fundamental and harmonic frequencies (described in greater detail below). The change in the IER (AIER) is estimated based on changes in those characteristics (e.g. magnitude and/or power).

[0039] The use of frequency domain processing for estimating AJER may provide a cost advantage over time domain signal processing. For example, the signal modulation caused by lung motion tends to have a relatively low frequency. As discussed in greater detail below, the time domain signal may be down-converted and filtered by narrow-band analog filters. Signals processed in these or other ways may be sampled at a relatively low rate (e.g. 128 Hz or even lower), thereby permitting the use of relatively inexpensive analog-to-digital converters and processors.

[0040] Figure 1A is a block diagram of a sensing system 14 for estimating ACER of a subject (not shown) according to a particular exemplary embodiment. Sensing system 14 may comprise a UWB radar system, a BCG sensing system, and/or another sensing system.

Sensing system 14 may estimate AIER of a human or animal subject. For brevity, but without loss of generality, the remainder of this description assumes that sensing system 14 is an IR- UWB sensing system used on a human subject, except where otherwise noted or where otherwise discernable from the context. Sensing system 14 may be used to detect a number of physiological characteristics associated with the subject. In addition to the change in the IER, these quantities include, without limitation, respiration rate (RR) and/or heart rate (HR).

[0041] In the illustrated embodiment, the operation of sensing system 14 is controlled by a suitably programmed controller 38. Controller 38 (and/or components thereof) may comprise hardware, software, firmware or any combination thereof. For example, controller 38 may be implemented on a programmed computer system comprising one or more processors, user input apparatus, displays and/or the like. Controller 38 may be implemented as an embedded controller with a suitable user interface comprising one or more processors, user input apparatus, displays and/or the like. Processors may comprise microprocessors, digital signal processors, graphics processors, field programmable gate arrays, and/or the like. Components of controller 38 may be combined or subdivided, and components of controller 38 may comprise sub-components shared with other components of controller 38. Components of controller 38, may be physically remote from one another. In some embodiments, controller 38 is not necessary and the operation of sensing system 14 may be controlled from a central monitoring station or the like which is remote from sensing system 14.

[0042] The example sensing system 14 of Figure 1A comprises a UWB radar system 30 connected to an antenna system 31. UWB radar system 30 comprises a UWB transmitter 3 OA, and a UWB receiver 3 OB. UWB radar system 30 generates UWB pulses (or impulses) that are transmitted by antenna system 31 into a space where the subject may be located. This space may be referred to as the "sensing volume" of sensing system 14, since UWB radar system 30 can be used to estimate the physiological characteristics of a person in the sensing volume. If a person is present in the sensing volume, then the pulses transmitted by UWB transmitter 30A are reflected at interfaces within the person's body (e.g. the surfaces of the lungs and heart) and the corresponding reflected pulses may be received at antenna 31 and detected by UWB receiver 30B.

[0043] In some embodiments, the UWB pulses used by UWB radar system 30 are in the C- band (e.g. in a range of 3 to 5 GHz). In some embodiments, the width of transmitted UWB pulses may be in a range of 1-100 ns, for example. UWB pulses are transmitted by UWB radar system 30 at a suitable pulse rate. This pulse rate may be set to a value low enough that the average power emitted by radar system 30 is low enough to satisfy applicable regulatory requirements. For example, the pulse repetition interval (PRI) is in the range of 0.5 to 1 in some embodiments. The time-averaged transmitted output power may be relatively small. For example, in some embodiments, the maximum effective isotropic radiation power (EIRP) may be -41.3dBm/MHz or less.

[0044] UWB radar system 30 generates a receiver output signal 32 based at least in part on reflected UWB pulses received by UWB radar receiver 30B. UWB radar system 30 may comprise signal conditioning components which are not expressly shown in Figure 1A for conditioning the signal received by UWB radar receiver 30B and thereby generating receiver output signal 32. Such signal conditioning components (e.g. circuitry, hardware and/or software) are well understood by those skilled in the art and are not described in detail here. By way of non-limiting example, such signal conditioning components may comprise amplifiers, filters, oscillators, mixers, gating circuits, analog to digital converters, sampling circuits and/or the like.

[0045] UWB radar system outputs a receiver output signal 32 which is passed to a signal processing system 33. Signal processing system 33 may comprise one or more digital signal processing (DSP) components configured with suitable DSP software for processing receiver output signal 32. In some embodiments, methods for estimating ΔΙΕΡν described herein may be implemented by signal processing system 33 and/or by signal processing system 33 under the control of controller 38 and may comprise using receiver output signal 32 to estimate AIER of a subject. Output signal 32 comprises measured sensor data from UWB radar system 30. In other embodiments, signal(s) comprising measured sensor data could come from other types of sensing systems (e.g. BCG systems and/or the like).

[0046] In the illustrated embodiment, signal processing system 33 provides an output 34 which is provided (either directly, or via controller 38) to optional input/output system 35. Output 34 may comprise information that is representative of the estimated physiological characteristics (including AIER) by sensing system 14, and may be provided, for example, in the form of a physical medium, change in state, and/or by any other means. In some embodiments, input/output system 35 may comprise a suitable user interface for providing output to users (e.g. visual output on a display and/or the like and/or audio output via speakers and/or the like) and/or for receiving input from users (e.g. via a keyboard, a pointing device, a touch screen interface and/or the like). User interfaces are well known and are not described further here. In some embodiments, for example where sensing system 14 is part of a monitoring system 10 (see Figure IB), input/output system 35 may comprise a

communication device (not expressly shown), such as a wireless communication device, which communicates a signal 24 comprising the estimated physiological characteristics (including AIER) to a remote monitoring station 12, by way of optional antenna 36.

[0047] Figure IB illustrates schematically how antenna system 31 of sensing system 14 located under a mattress 17 can emit UWB pulses 40, some of which pass into a subject S and are reflected at various surfaces, including surfaces of the subject's lungs L and heart H, to form reflected pulses 42 that are in turn detected by antenna system 31. By using UWB radar system 30 in a back-scattering remote sensing mode, sensing system 14 can estimate physiological characteristics (e.g. heart and respiration rates) of a subject S without electrodes or other devices being attached to the body of subject S.

[0048] Antenna system 31 may comprise an array of transmit antennas (not expressly shown) and an array of receive antennas (not expressly shown). The transmit and receive antennas may be distributed over an area broad enough to be able to transmit UWB pulses 40 and detect reflected UWB pulses 42 from subject S in any reasonable position and posture on mattress 16. The transmit antennas may be low-gain antennas (relative to the receive antennas). The use of low-gain transmit antennas permits UWB signals having higher average amplitudes to be applied to the antenna without causing EIRP to exceed thresholds that maybe specified by applicable regulations. Further, low-gain antennas generally have broad radiation patterns the radiation is distributed into a broad angular sensing volume. The receive antennas of antenna system 31 may comprise higher-gain antennas (than the transmit antennas) to provide better signal-to-noise ratios (SNR) for received signals ascertained from reflected pulses 42.

[0049] While not necessary, sensing system 14 may be part of a monitoring system (e.g. a patient monitoring system). Some such monitoring systems are described in US Patent Application No. 14/247070, referenced above and previously incorporated by reference.

[0050] Figure 2A is a block diagram representation of a method 110 for estimating ATKR of a subject according to a particular embodiment. Method 110 may be performed by sensing system 14 (e.g. by digital signal processing system 33 or digital signal processing system 33 in combination with controller 28). In some embodiments, method 110 may be performed in part by UWB radar system 30 under the control of signal processing system 33 and/or controller 28. Alternatively, or in addition, method 110 may be performed in part by a BCG sensing system, substantially as described herein.

[0051] Method 110 starts in block 112 which comprises obtaining frequency domain data 1 16. In the illustrated embodiment, block 112 comprises processing time domain

measurement data 114 output from a UWB radar system (e.g. UWB radar system 30 of sensing system 14 of Figure 1A) to obtain a frequency domain representation of time domain measurement data 114 which may be output as frequency domain data 116. When UWB transmitter 30A transmits pulses (e.g. pulses 40 shown in Figure IB) with a period T, these pulses are reflected (e.g. reflected pulses 42) from the various surfaces of the body of the subject and have a time varying delay which varies with movement of the various surfaces of the body and may be modelled by the equation: r(t) = a _h p{t - nT - _Tft (t)) + a _lV (t - nT - r,(t)) + a _s p(t - nT - T _s (t)) n=-∞ where the first term in the sum is representative of reflections from the heart of the subject, the second term in the sum is representative of reflections from the lungs of the subject and the third term in the sum is representative of reflections from skin of the subject. In equation (1): p(. ) represents a particular pulse, Tis the pulse repetition period, e _s are amplitudes of the reflected pulses coming from heart, lung and skin interface while τ _¾ , Tj and T _s is the time delay of the reflected pulses coming from heart, lungs and skin.

[0052] In some embodiments, UWB radar system 30 transmits a periodic series of brief, low- power pulses (e.g. pulses of a 1-100 nanoseconds transmitted with a period of 1-4μ8). When those pulses are reflected by a surface that is in approximately periodic (e.g. oscillatory) motion, their time of arrival is altered relative to the average delay, causing a delay modulation of the reflected signal. For example, when the lung's surface moves towards the IR-UWB radar system 30, the reflected pulses 42 are received sooner, while when the lung's surface moves away, the reflected pulses 42 are delayed. The resulting delay-modulated time domain signal may then be analyzed to identify the frequency of the oscillatory motion inducing the delay (e.g. the respiration rate associated with the movement of the lung).

[0053] Time domain measurement data 114 may comprise a version of the radar-received signal received by UWB radar receiver 30B which may be down-converted to baseband. As is known in the art, such down-conversion may be performed by suitable components (not expressly shown) of UWB radar receiver 30B. In some embodiments, after down-conversion, time domain measurement data 114 may be pre-filtered by analog filters (not expressly shown) prior to being received in block 112. Such analog pre-filters may be a part of UWB radar receiver 30B and may comprise band-pass filters which attempt to separate the fundamental frequency of one physiological characteristic (e.g. the respiration rate) from the fundamental frequencies of other physiological characteristics (e.g. the heart rate). Some such pass-band filtering processes are described in greater detail in US Patent Application No. 14/247070, referenced above and previously incorporated by reference.

[0054] Figure 2B schematically illustrates a pair of band-pass filters 126, 128 which are exemplary of the type of filters which may be used to filter respiration data from heart rate data. In some embodiments, the signal received by block 112 may be the band-pass filtered signal 130 which is passed by the respiration filter 126. Accordingly, in some embodiments where the physiological characteristic of interest is the respiration rate, time domain measurement data 114 received in block 112 of Figure 2A may comprise this band-pass filtered signal 130 shown in Figure 2B (or a digitally sampled version thereof). Respiration filter 126 may be designed to have a pass-band which captures the respiration rate fundamental and one or more respiration rate harmonics for a suitable range of respiration rates. [0055] In some embodiments, time domain measurement data 114 received from UWB radar system 30 may be sampled and digitized prior to being received in block 112 (e.g. by suitable digital sampling components (not expressly shown) which may form part of UWB receiver 30B). Such digital sampling components are well known in the art. Unlike some prior art techniques which involve sampling rates of tens of GHz, the sampling of time domain measurement data 114 may be performed at a relatively low rate. In one particular embodiment, the sampling rate of time domain measurement data 114 is around 32Hz. In some embodiments, this sampling rate is less than or equal to 128Hz. In some embodiments, this sampling rate is less than or equal 256Hz. In some embodiments, this sampling rate may be higher if desired (for example, to use particular hardware). In some embodiments, time domain measurement data 114 received from UWB radar system 30 may still be in the analog domain and may be sampled and digitized as a part of block 112. The block 112 sampling rates may be comparable to those where the sampling and digitizing is performed by UWB receiver 30B.

[0056] In the illustrated embodiment, block 112 comprises converting time domain measurement data 114 into frequency domain data 116. Frequency domain data 116 may comprise an estimate of the power spectral density (PSD), amplitude spectral density or other suitable spectral density estimation of the reflected pulses received by UWB receiver 30B (and the corresponding time domain measurement data 114). In general, any suitable technique may be used to convert or transform time domain measurement data 114 to the frequency domain to obtain frequency domain data 116. The particular technique used to generate frequency domain data 116 may depend on the nature of time domain measurement data 114. The frequency domain data 116 generated from a set of time domain measurement data 114 may be referred to herein as a frequency domain representation of that time domain measurement data 114.

[0057] A challenge with measuring the movement of subjects' organs or other surfaces using UWB radar system 30 is that time domain measurement data 114 may be sensitive to other movements of the subject. Persons skilled in the art will appreciate that the delay modulation induced by a subject simply moving a limb may be orders of magnitude greater than the delay modulation caused by the movement of the subject's lungs within the subject's body, for example. Specific embodiments and implementations of particular features of method 110 directed to improving sensing system 14's accuracy and/or reliability in at least some circumstances are now described in more detail.

[0058] Figure 5A is a block diagram representation of a method 200 for obtaining frequency domain data 116 according to a particular embodiment. In some embodiments, method 200 of Figure 5 A may be used to implement block 112 of method 110 (Figure 2A). Method 200 may be performed by UWB sensing system 14 (e.g. by digital signal processing system 33 or digital signal processing system 33 in combination with controller 28). In some embodiments, method 200 may be performed in part by UWB radar system 30 under the control of signal processing system 33 and/or controller 28.

[0059] Method 200 commences in block 202 which involves receiving time domain data 114 from a UWB radar system (e.g. UWB radar system 30 (Figure 1A)). As discussed above, time domain data 114 may comprise a version of the radar-received signal received by UWB receiver 30B and down-converted to baseband. As also discussed above, time domain data 114 may be pre-filtered by analog filters to separate the respiration rate fundamental from other data (such as the heart rate fundamental) and to block DC.

[0060] Block 202 of the illustrated embodiment comprises sampling time domain data 114 at a first sampling rate ri. This block 202 first sampling rate ri may be relatively low as compared to prior art sampling rates. In some embodiments, this first sampling rate is 32Hz. In some embodiments, this first sampling rate is less than or equal to 128FIz. In some embodiments, this first sampling rate is less than or equal 256FIz.

[0061] In block 204, an optional noise suppression process is applied to the data sampled in block 202. For example, block 204 may apply a DC removal filter and/or a Hamming window filter to time domain measurement data 114 prior to transforming time domain measurement data 114 to the frequency domain.

[0062] In some embodiments, block 204 comprises applying a down-sampling or decimation process to the data sampled in block 202 to generate down-sampled (decimated) time domain data. This block 204 down-sampling process may reduce noise which may be present in the sampled data output from block 202. The optional block 204 down-sampling process may involve taking the average of every j samples of the block 202 data and generating a single sample in the block 204 down-sampled data set for every j samples of the block 202 data. In some embodiments, j ^* =8, although j could be provided with other values. The block 204 down-sampling effectively reduces the block 202 sampling rate to Thus, in example embodiments where Hz and y-8, Γ2 is 16 Hz. In some embodiments, the optional noise suppression functionality of block 204 may be provided by different types of operations.

[0063] In some embodiments, block 204 may attempt to isolate acute motion events by examining time domain measurement data 204, the data sampled in block 202, and/or the decimated time domain data (collectively and individually the "sample data") for irregularly large values. Block 204 may, for example, perform DC removal filtration on the sample data prior to checking for irregularly large values to avoid incorrectly identifying typical values as irregularly large due to subject motion. Block 204 may provide a moving average filter for this purpose. For example, block 204 may determine an average amplitude of the sample data over a relatively long period {e.g. over 10, 20, 30, 45, 60, 90, or 120 seconds or more) and remove that average from the sample data (e.g. by subtraction) to filter out DC bias.

Additionally, or in the alternative, block 204 may determine the absolute value of the sample data to include large negative values in its examination.

[0064] Once block 204 has suitably adjusted the sample data (if any adjustment is performed), block 204 may examine the sample data, as adjusted, for irregularly large values and omit those values so that they are not included in the conversion or transformation from the time domain into the frequency domain. For example, block 204 may compare amplitudes of the sample data to a threshold and, if those amplitudes are greater than the threshold, those amplitudes may be omitted. A suitable threshold may, for example, be configurable {e.g. based on experimental observation or heuristics and/or the like), determined based on an average {e.g. the threshold may be set to 200% of the average value of amplitudes over the last 60 seconds), or otherwise determined {e.g. the threshold may be set to 400% of the bottom quartile of the amplitudes time domain measurement data 114 over the last 120 seconds). In some embodiments, the omission of amplitudes greater than the threshold is performed elsewhere {e.g. in any of blocks 206, 208 and/or 210).

[0065] Other techniques may additionally or alternatively be used to remove spurious data from time domain measurement data 114. By way of non-limiting example, short-term estimates of signal variance may be determined and smoothed (using any suitable technique) and data for which the smoothed variance is above a threshold may be excised. As another example, a nonlinear operation may be performed on the data (e.g. squaring) to increase contrast, before smoothing and thresholding. In some embodiments, a second threshold may be applied to the number or percentage of time-domain samples that exceed the first amplitude threshold. For example, if more than 5% of the sample magnitudes (after subtracting the DC trend) exceed a certain threshold, the spectra associated with that time interval may be considered "motion corrupted" and not used in further processing.

[0066] By way of non-limiting example, in some embodiments, the block 204 noise suppression functionality may be provided by a moving average filter, such as described above. Such a moving average filter may not involve a reduction in the effective sampling rate of the sample data (although, in some embodiments, a moving average filter may induce a delay in providing at least some filtered data to block 206). In other embodiments, block 204 may comprise applying other suitable averaging or down-sampling filters to the data sampled in block 202 to suppress noise. In some embodiments, block 204 is not necessary.

[0067] Method 200 then proceeds to block 206 which involves parsing the sampled time domain data (output from block 202 or from optional block 204) into temporal windows. In some embodiments, block 202 may receive time domain measurement data 114 in the form of multiple temporal windows. Alternatively, or in addition, method 200 may parse time domain measurement data 114 into temporal windows at block 206 (e.g. block 206 may subdivide temporal windows or determine new temporal windows which are not necessary bounded by received temporal windows, if any).

[0068] Figure 5B is a schematic depiction of a method 220 for parsing time domain measurement data into temporal windows which may be used in block 206 according to a particular embodiment. Method 220 of Figure 5B involves parsing the sampled time domain data into temporal windows of duration tj, where each successive window is offset from the preceding window by an offset t ₀ , where For example, in accordance with method 220, a first block 206 temporal window Wi may include samples of the time domain data output from block 202 or from optional block 204) during the period between t=0 and t= , a second temporal window W ₂ includes samples of the time domain data during the period between t=to and a third temporal window W3 includes samples of the time domain data during the period between t=2to and t=tj+2to, a fourth temporal window W4 includes samples of the time domain data during the period t=3to and t=ti+3to and so on.

[0069] With this method 220 scheme, each successive temporal window shares samples (of a duration tj-to) with its preceding temporal window. In some embodiments, the ratio of a duration of shared samples between successive windows to a duration of each window (i.e. ) is in a range of 50% to 90%. In one particular embodiment, this ratio is 75%. In some embodiments, the duration of each window tj is selected to be in a range of 5-60 seconds. In one particular embodiment, tj is set to be seconds. In some embodiments, the offset to is selected to be in a range of 0.5-30 seconds. In one particular embodiment, ¾ is set to be t ₀ =4 seconds. For ease of explanation, it will be assumed in the remainder of this description (without loss of generality) that the block 206 windowing is configured such that seconds, to=4 seconds and the ratio is 75%. Further, it will be assumed without loss of generality that the sampling rate r _/ =128Hz and that above-described block 204 down- sampling occurs with j=8. With these assumptions, each temporal window W _1} W ₂ , W3, . . . includes 256 samples.

[0070] Returning to Figure 5A, after parsing the time domain data into finite duration windows in block 206, method 200 proceeds to block 208 which comprises optional pre- transform processing. In some embodiments, the optional block 208 pre-transform processing comprises DC removal and Hamming windowing. Hamming windowing may increase the relative amplitude of the main lobe of spectral components and decrease relative amplitude of the side lobes, which may in turn reduce interference between spectral components. In some embodiments, other types of windowing (e.g. Harming windowing, Kaiser windowing and/or the like) and/or other types of pre-transform processing may be performed in block 206. This processing may be applied independently to the time domain data in each of the block 206 windows.

[0071] Method 200 then proceeds to block 210 which involves transforming the time domain data output from block 208 to the frequency domain. In general, any of the techniques described above in relation to blocks 112, 202, and/or 204 may be applied independently to each block 206 temporal window Wi, W ₂ , W3, ... of time domain data. In some

embodiments, the block 210 transformation involves determining an estimated power spectral density (PSD), amplitude spectral density or using some other suitable spectral density estimation technique for each block 206 temporal window Wi, W ₂ , W3, ... of time domain data. For example, in one particular embodiment, block 210 comprises applying a DFT to each block 206 temporal window Wi, W ₂ , W3, ... of time domain data and squaring the magnitude of the result to obtain a PSD estimate for each block 206 temporal window.

[0072] As is the case with block 112 discussed above, any suitable technique may be used to compute the DFT, such as, by way of non-limiting example, any suitable FFT technique. In other embodiments, other suitable spectral transform techniques or other spectral density estimation techniques could be used in block 210 in the place of the magnitude-squared DFT. Such spectral transform techniques include, by way of non-limiting example, any suitable spectral estimation technique (e.g. MUSIC, ESPRIT or the like), any suitable time-frequency transform technique, any suitable transform involving a weighted sum of orthogonal and/or periodic functions, and/or the like. In some embodiments, it is not necessary to square the amplitude of the transformed data. For example, in some embodiments wherein time domain measurement data 114 is obtained from a BCG sensor the amplitude of the transformed data need not be squared.

[0073] As described above, each temporal window of time-domain data may be transformed into a frequency domain representation, thus providing a plurality of frequency domain representations wherein each frequency domain representation is associated with a corresponding temporal window. In general, unless the context provides otherwise, any of the techniques applied to frequency-domain data in block 112 of method 110 and/or blocks 210, 212 and 214 of method 200 may be applied independently to each of these frequency domain representations.

[0074] In some embodiments, block 210 is the last step of method 200 and results in frequency domain data 116. In the illustrated embodiment, method 200 involves a number of optional procedures shown in blocks 212 and 214. The procedures of blocks 212 and 214 are not required. In some embodiments, the procedures of block 212 may be performed without the procedures of block 214 or the procedures of block 214 may be performed without the procedures of block 212. In other embodiments, a more finely sampled version of the spectrum could be obtained by other techniques. For example, such other techniques may involve zero-padding the time domain data from the block 206 temporal windows Wi, W ₂ , W3, ... and/or the like. In other embodiments, a higher resolution version of the spectrum could be obtained by other techniques. For example, such other techniques may involve use of longer temporal windows Wi, W ₂ , W3, ... (i.e. windows with a larger parameter ti).

[0075] Optional block 212 involves averaging the block 210 frequency domain data for a plurality of successive temporal windows Wi, W ₂ , W3, ... (i.e. windows Wj, W ₂ , W3, . . . obtained in block 206). In one particular embodiment, this block 212 averaging is performed for each pair of successive temporal windows W _ls W ₂ , W3, This particular embodiment of block 212 is depicted schematically in the method 222 of Figure 5B. In the Figure 5B illustration, the vertical dashed line represents the block 210 transformation of time domain data to frequency domain data. Each block 206 temporal window results in corresponding estimated power spectral density data PSD(Wi), PSD(W ₂ ), PSD(W ₃ ) in the frequency domain. In block 212, successive pairs of frequency domain data are averaged - i.e. PSD(Wi) and PSD(W ₂ ) are averaged to form PSD(Wi,W ₂ ), PSD(W ₂ ) and PSD(W ₃ ) are averaged to form PSD(W ₂ ,W ₃ ), PSD(W ₃ ) and PSD(W ₄ ) are averaged to form PSD(W ₃ ,W ₄ ) and so on. In other embodiments, the spectral densities (frequency domain data) of different numbers of temporal windows could be averaged in block 212.

[0076] Method 220 then proceeds to optional block 214 which involves interpolating and/or smoothing the averaged spectra from the block 212 averaging process. Block 214 may provide an approximation of a more finely sampled version of the spectrum. Any suitable interpolation and/or smoothing technique could be used in block 214. In one particular embodiment, block 214 comprises applying a cubic spline interpolation technique. In some embodiments, forward and backward recursive filtering may be combined to produce a smoothed version of the spectrum sequence. Where the block 212 averaging technique involve the averaging of pairs of successive temporal windows to generate frequency domain data PSD(Wi,W ₂ ), PSD(W ₂ ,W ₃ ), PSD(W ₃ ,W ₄ ) ... (as in the case in the illustrated embodiment of Figure 5B), the block 214 interpolation may be applied to these frequency domain data to obtain interpolated and averaged frequency domain data PSD(W1, W2), PSD(W2, W3), PSD(W3, W4) .... This interpolated and averaged frequency domain data PSD(W1, W2), P!TD(W2, W3), Pl¾(W3, W4) .... may be output as frequency domain data 116 shown in Figure 5A - i.e. the result of method 200 and block 112 (Figure 2A). For brevity but without loss of generality, the remainder of this description, except where otherwise noted or where otherwise discernable from the context, assumes that frequency domain data 116 has the form of interpolated and averaged frequency domain data

P D(W1, W2), PSb(W2,W3), PSD(W3, W4) ... output from block 214.

[0077] The PSD of a periodic signal has peak values at or near various frequencies. Of particular interest in the present disclosure are the fundamental frequency (sometimes referred to simply as the "fundamental" and/or denoted f), which is the inverse of the period of the signal, and the integer multiples of /(sometimes referred to as the "harmonic frequencies", or simply "harmonics", of f). The harmonics of a particular fundamental may be identified as the "first harmonic" (e.g. RRi in Figure 3), "second harmonic" (e.g. RR ₂ in Figure 3), and so on.

[0078] Figure 3 shows an example of the type of frequency domain data 116 that could be obtained from block 112 according to some embodiments. The data shown in the particular example of Figure 3 may be obtained in block 112 by application of a DFT to a 20 second window of digitally sampled time domain data and then squaring the magnitude of the DFT to obtain the PSD estimate, although, as discussed above, data similar to that shown in Figure 3 could be obtained using other techniques. The unit shown on the horizontal axis of the Figure 3 plot is frequency. As is known in the art, the DFT will commonly produce "bins" of data, which correspond generally to frequency values. Those skilled the art will appreciate that there is a relationship between DFT bins and the frequency of the digitally sampled time domain data which depends on the time domain sampling rate (which, in this example, is 64 FIz) and the length of the window over which the DFT is obtained (which, in this example, is 20 seconds). For the sake of intelligibility, frequency is shown along the horizontal axis in Figure 3, but those skilled in the art will appreciate that the systems and methods described herein may in some circumstances be used without necessarily determining the frequency values associated with particular DFT bins.

[0079] In the illustrated data shown in Figure 3, the respiration rate in breaths per minute (bpm) is approximately 8, indicated by the position of the highest peak (RRo) along the frequency axis. RRo is the fundamental frequency of the signal associated with respiration rate. For convenience, Figure 3 shows the fundamental respiration rate RRo as a dashed vertical line with a closed arrowhead and the various harmonics RRi, RR ₂ , RR ₃ as dashed vertical lines with open arrowheads. The respiration rate harmonics RRi, RR ₂ , RR3 shown in Figure 3 may be approximately integer multiples of the respiration fundamental RRo.

Frequencies RRo, RRi, RR ₂ , and RR ₃ have magnitudes P ₀ , Pi, P ₂ , and P ₃ , respectively.

[0080] In some embodiments, obtaining frequency domain data 116 in block 112 involves merely receiving frequency domain data 116 (rather than processing time domain measurement data 114 to obtain frequency domain data 116). In such embodiments, frequency domain data 116 may be received as a part of block 112 from some other source (e.g. from some other measurement or sensing system, from some other processor and/or the like) which may compute frequency domain data 116 and provide same to block 112. [0081] Once frequency domain data 116 is obtained in block 112, method 110 (Figure 2A) proceeds to block 118 which involves using frequency domain data 116 to discern one or more characteristics of the frequency domain data 116 at respiration fundamental RRo and the one or more respiration harmonics RRi, RR ₂ , RR ₃ , .... For example, block 118 may discern the magnitude, frequency, power, and/or other values associated with the fundamental and/or harmonic(s). In general, a characteristic of frequency domain data 116 at the fundamental frequency may be referred to herein as a "fundamental characteristic". Particular types of characteristics may be referred to similarly; for example, the magnitude of frequency domain data 116 at the fundamental frequency may be referred to herein as the "fundamental magnitude". Similarly, a characteristic of frequency domain data 116 at a harmonic frequency may be referred to herein as a "harmonic characteristic". A particular characteristic of frequency domain data 116 at a particular harmonic may be identified as, for example, the "first harmonic magnitude", the "second harmonic power", "the third harmonic

neighbourhood power", or the like.

[0082] The discerned characteristics of the fundamental and/or harmonic(s) are the result or output of block 118, and may be referred to herein as discerned characteristics 120. Unless the context dictates otherwise, the discerned characteristics 120 pertain to the respiration fundamental and one or more of its associated respiration harmonics (although, optionally, discerned characteristics 120 may additionally comprise characteristics of the heart rate fundamental, its harmonics, and/or other fundamentals and their harmonics).

[0083] Discerning fundamental and harmonic characteristics may involve detecting "peaks" in the frequency domain data 116 associated with a fundamental and/or its harmonics. Some approaches to identifying such peaks, including harmonic path detection of peaks, are discussed in greater detail in US Patent Application No. 14/247070, referenced above and previously incorporated by reference. Alternatively, or in addition, the fundamental and/or harmonics may be identified using any other technique (e.g. by estimation of the respiration rate via a capnograph, Respiband, pneumotachograph or the like), and the peaks may be identified based on the identified fundamental and/or harmonic frequencies.

[0084] In the case of the illustrative Figure 3 frequency domain data 116, it can be seen that frequency domain data 116 exhibits frequency domain peaks at frequency locations which correspond to the respiration fundamental RRo and at least a number of the respiration harmonics RRi, RR ₂ , RR3. Accordingly, in some embodiments, block 118 comprises discerning ascertaining local maxima (peaks) in frequency domain data 116. In some embodiments, such peak detection may comprise application of a suitable identification threshold to select peaks of interest (e.g. peaks which may be considered to be located at candidates for the respiration rate fundamental and/or harmonics of the respiration rate may be discerned to have local maxima that are greater than the identification threshold). An exemplary identification threshold 132 that may be used in the block 118 thresholding process is shown as a horizontal line in Figure 3.

[0085] In some embodiments and/or for some particular data, the block 118 thresholding process may be sufficient to discern a plurality of peaks corresponding to a fundamental and/or its harmonics (e.g. by suitable selection of a threshold). In the example case illustrated by Figure 3, identification threshold 132 assists in discerning the respiration fundamental, but not any of its harmonics. In some embodiments, it may be desirable to set the identification threshold of the block 118 thresholding process to be somewhat lower than the exemplary identification threshold 132 shown in Figure 3 (e.g. to avoid accidental exclusion of suitable peaks or to discern a greater number of candidate peaks).

[0086] An example of a relatively low identification threshold 134 that may be used in the block 118 thresholding process is shown as a horizontal line in Figure 3. Where threshold 134 is applied to frequency domain data 116 as a part of block 118, an additional peak (of the Figure 3 exemplary data) satisfies the block 118 thresholding criteria. Consequently, in some embodiments, it may be desirable to perform one or more additional or alternative evaluation processes (e.g. to evaluate one or more additional or alternative criteria) as a part of block 118 before outputting the discerned characteristics 120. Such additional or alternative evaluation processes may help to more accurately discern a plurality of elements from the respiration rate harmonic set. Some such additional or alternative processes are described in greater detail in US Patent Application No. 14/247070, referenced above and previously incorporated by reference.

[0087] In some embodiments, the block 118 identification threshold is set as some percentile of the cumulative distribution function (CDF) of the amplitudes of frequency domain data 116. For example, the block 118 identification threshold may be set to the 75 ^th percentile of the amplitudes of frequency domain data 116. Other suitable CDF percentile thresholds may be used and such percentile thresholds may be configurable (e.g. user-configurable). In other embodiments, other criteria (such as experimental data) may be used to configure the identification threshold used to detect peaks in block 118. Alternatively, or in addition, a identification threshold may be set relative to the largest amplitude in frequency domain data 116 (e.g. 10%, 25%, 50%, 75%, or 90% of the largest amplitude) and/or relative to an average value of the signal (e.g. 100%, 125%, 150%, 175%, 200%, 250%, or 300% of an average value). A identification threshold may be a constant value, or it may be a function that varies with the frequency (or DFT bin) domain (i.e. different frequencies or DFT bins within frequency domain data 116 may be compared against different values during the thresholding process).

[0088] In embodiments which use average values to determine a identification threshold, such average values may be determined in any of several ways. In some embodiments, the average value may be determined by taking an average of some or all of the magnitudes of frequency domain data 116. For example, the average value associated with a particular frequency or DFT bin (i.e. the value based on which the threshold for that bin will be determined) may be an average of the magnitudes of nearby frequencies (or DFT bins). For example, identification threshold 136 is an example of a identification threshold based on an average value that, at each frequency, is based on the average of all magnitudes in a neighbourhood of 20bpm. For instance, the value at bpm=8 in the exemplary identification threshold 136 of Figure 3 is based on the average of the magnitudes in the interval [0,18], and the value at bpm=20 is based on the average of the magnitudes in the interval [10,30].

Identification threshold 136 may be used for the thresholding processes described herein, in some embodiments.

[0089] In some embodiments, the average values used to determine identification threshold 136 are determined based on multiple sets of frequency domain data 116. For example, the average value may be an average of the sets of frequency domain data 116 associated with a consecutive series of windows of time domain measurement data 114 (e.g. based on the current window and one or more windows preceding it). For instance, an average value may be determined by taking the average of all of the magnitudes of some or all of the frequency domain samples in the previous three windows of frequency domain data 116. Such an average can be referred to as a moving average.

[0090] In some embodiments, the fundamental frequency is discerned by identifying the frequency in the frequency domain data 116 with the greatest magnitude. [0091] In some embodiments, one or more harmonic frequencies are identified by determining integer multiples of the fundamental frequency (regardless of how the fundamental frequency is identified). For example, for a fundamental frequency the first harmonic frequency (e.g. RRi in Figure 3) may be determined to be 2f or the highest value in a range of 2/ ± Δ, the second harmonic frequency may be determined to be 3/ or the highest value in a range of 2/ + Δ, and so on, where Δ is a configurable parameter. In some embodiments, this concept (i.e. that harmonics are located at or near integer multiples of fundamental frequency f) may be used as a verification test for harmonic peaks acquired using other techniques.

[0092] The above techniques for identifying the fundamental and/or harmonic frequencies and/or their associated peaks are provided as non-limiting examples only. As noted above, the fundamental and/or harmonic frequencies may be identified using any other technique. For example, the fundamental and/or harmonic frequencies may be identified by estimation of the respiration rate via a capnograph, Respiband, pneumotachograph or the like. Such techniques may not necessarily require identification of peaks associated with the fundamental and/or harmonic frequencies.

[0093] In some non-limiting embodiments, discerned characteristics 120 comprises only the fundamental and the first harmonic frequency of the physiological characteristic(s) of interest. For example, sensing system 14 may discern only the fundamental and first harmonic frequencies of the respiration rate.

[0094] Once discerned characteristics 120 are determined in block 118, method 110 proceeds to block 122 which involves using discerned characteristics 120 determine a AIER estimate 124. Block 122 may comprise determining AIER estimate 124 based on a change (or multiple changes) in the discerned characteristics 120. For example, the exemplary Figure 3 data has a magnitude P ₀ at fundamental frequency RRo and magnitudes Pi, P ₂ and P ₃ at harmonics RRi, RR ₂ and RR ₃ , respectively. The magnitudes P ₀ , Pi, P ₂ , and P ₃ may be referred to as magnitudes of their respective frequencies. For example, assuming that at least the peaks located at RRo and RRi were discerned by block 118, block 122 may determine changes in the discerned characteristics associated with those peaks (e.g. P ₀ and Pi) relative to one or more previously-analyzed sets of frequency domain data 116 (e.g. sets of frequency domain data 116 associated with one or more previous temporal windows of time domain measurement data 114). Block 122 may then determine AIER estimate 124 based on those relative changes over time.

[0095] AIER estimate 124 produced by block 122 may comprise an estimate of a cumulative change in IER (e.g. an estimate of the total change in IER since sensing system 14 began monitoring the subject), and/or AIER estimate 124 may comprise a discrete change in IER (e.g. an estimate of the change in IER since the previous AIER estimate 124 was provided). Without loss of generality, AIER estimate 124 will be assumed to be cumulative in the following disclosure, unless the context provides otherwise. In some embodiments, sensing system 14 may estimate IER, for example by obtaining an initial measurement of IER (e.g. via measurement by a capnograph, Respiband™, and/or pneumatotachograph) and estimating IER based on that initial measurement and subsequent ER estimate(s) 124.

[0096] The inventors have observed, through analysis, that a positive change in IER (i.e. increasing IER) tends to correlate with a relative growth in the peak of the respiration fundamental and/or a relative shrinkage of the peaks of the respiration harmonics. Without wishing to be bound by any particular explanation, this phenomenon may be explained as follows. When the bronchial passages become more inflamed, muscles work harder to pull the air in (this is why a child becomes exhausted), but the expiration, which is normally a passive process, takes a longer time, making the IER smaller. As IER decreases, slope of the lung or skin displacement function during the inspiration period increases, when the function is plotted against normalized time (normalized by the respiration period). Higher frequency Fourier components may describe this higher slope. If the overall norm (e.g., energy) of the displacement function is conserved, then the power in the PSD may shift from the fundamental out to the harmonic components.

[0097] Note that most published analysis of Fourier components of the IR-UWB signal directed towards estimation of respiration and/or heart rate assume perfect sinusoidal displacement, which has no higher frequency components in the PSD of the displacement. The IR-UWB signal, which is delay-modulated by the displacement function, similarly reflects this shift of power from fundamental to harmonic components. The BCG signal is approximately proportional to displacement, by definition, and so its spectrum also exhibits this shift of power from fundamental to harmonic components as IER changes. Accordingly, in some embodiments, discerned characteristics 120 may comprise the power of the PSD (and/or another frequency domain representation of the signal) in a neighbourhood of the fundamental frequency and/or the harmonic frequencies. This is discussed in greater detail below, along with the notion of estimating the power of the signal.

[0098] Note also that a signal (e.g. as represented by frequency domain data 116) may have multiple fundamentals, e.g. a respiration fundamental (corresponding to the respiration rate of the subject), a heart rate fundamental (corresponding to the heart rate of the subject), and so on. Multiple sets of fundamentals and/or harmonics may intermodulate, producing additional peaks which are sometimes referred to as "intermods". In some circumstances, the characteristics (e.g. magnitude, power, etc.) of intermods are proportional to the

characteristics of the fundamentals and/or harmonics which intermodulated to give rise to the intermods. For example, an intermod generated by the intermodulation between the heart fundamental and a respiration harmonic will tend to double in magnitude if the respiration harmonic magnitude doubles. Accordingly, since characteristics of intermods may change with the characteristics of respiration fundamentals and/or harmonics, in some embodiments intermods are also (or alternatively) analyzed by block 122 in its determination of ATER estimate 124.

[0099] Figure 4A shows plots of example experimental data in the frequency domain derived from an example IR-UWB signal (e.g. as shown in Figure 1A). In each of Figures 4A and 4B, example experimental data comprises data derived from simulations modelling simulated lung interfaces and JR-LTWB and/or BCG signals. Figure 4A shows a fundamental plot 150a showing the magnitude of frequency domain data 116 at the fundamental frequency and a harmonic plot 150b showing the magnitude of frequency domain data 116 at the first harmonic. These magnitudes are identified in Figure 3 as Po (for the fundamental) and Pi (for the first harmonic), respectively, although the exemplary Figure 4A does not necessarily use the same data as the exemplary Figure 3. In each plot 150 the horizontal axis shows IER and the vertical axis shows a measure of the magnitude. In the case of fundamental plot 150a, the measure used is simply the fundamental magnitude (for example, as determined according to the PSD estimate). In the case of harmonic plot 150b, the measure used is the normalized first harmonic magnitude, i.e. the magnitude of the first harmonic as a proportion of the magnitude of the fundamental, P ₀ .

[0100] The values modelled on the vertical axes also vary with displacement of the surface of reflection - e.g., in the case of the respiratory fundamental and harmonics, the surface of the subject's lung. Data 152a, 154a, 156a, and 158a correspond to a characteristic (in this example, the peak P _o ) of the fundamental in experiments wherein the subject's lung was displaced by ±2mm (i.e. 4mm of total displacement between the nearest and furthest positions of the lung's surface relative to sensing system 14), ±4mm, ±6mm, and ±8mm, respectively. Data 152a, 154a, 156a, and 158a show that, as IER increases, so does the magnitude of frequency domain data 116 at the fundamental frequency, although the rate of increase varies depending on the displacement of the subject's lung.

[0101] Data 152b, 154b, 156b, and 158b correspond to the ratio Pi P ₀ in the same experiments as described above. That is, data 152b, 154b, 156b, and 158b show the first harmonic magnitude as a proportion of the fundamental magnitude in frequency domain data 116. As shown, for IERs between 0.15 and 0.5, the ratio Pi/Po tends to grow with increasing displacement, thereby obscuring the relationship with IER.

[0102] Similarly, Figure 4B shows plots of example experimental data in the frequency domain derived from an example BCG signal. Namely, Figure 4B shows a fundamental plot 160a showing the magnitude of frequency domain data 116 at the fundamental frequency and a harmonic plot 160b showing the magnitude of frequency domain data 116 at the first harmonic. As in Figure 4A, each plot 160 shows IER along the horizontal axis and a measure of magnitude along the vertical axis. In the case of fundamental plot 160a, the measure used is the fundamental magnitude (e.g. as determined according to the PSD estimate) normalized to its smallest value (i.e., in this example, the leftmost value is set to 1 on the vertical axis). As in Figure 4A, the measure used in harmonic plot 160b is the normalized first harmonic magnitude, i.e. the magnitude of the first harmonic as a proportion of the magnitude of the fundamental.

[0103] As in Figure 4A, data 162a correspond to a characteristic (e.g. the peak P ₀ , normalized as described above) of the fundamental, and data 162b correspond to a ratio of a fundamental characteristic and a first harmonic characteristic (e.g. the ratio Pi/Po). Unlike Figure 4A, additional data are not provided in plots 160a and 160b for various chest displacements. This is because the BCG signal varies linearly with chest and/or lung displacement, unlike the IR- UWB signal (which varies non-linearly with chest and/or lung displacement). Accordingly, the frequency domain representation of the IR-UWB signal also varies linearly with chest and/or lung displacement. By normalizing the magnitudes shown in plots 160a and 160b, this linear variation is cancelled out.

[0104] In the example data shown in Figures 4A and 4B, the slope of the plot 150b and 160b data (e.g., the slope of the ratio Pi/Po) remains relatively constant under certain

circumstances. For example, in plot 150b, the slope is relatively constant for a given lung displacement for IERs in the interval 0.15 to 0.5. As another example, in plot 160b, the slope is relatively constant for IERs in the interval 0.2 to 0.5 (and the slope outside of this interval shows a similar downward trend) Accordingly, a change in the ratio of a characteristic of the frequency domain data 116 at the fundamental frequency and a characteristic of the frequency domain data 116 at a harmonic frequency (e.g. the ratio Pi/Po) may be used to estimate a change in IER (even in circumstances where IER may not be estimated directly). In some embodiments, harmonic characteristics other than a first harmonic characteristics (e.g. characteristics associated with RR ₂ , RR ₃ , ...) may be used in addition to, or instead of, the first harmonic characteristic, as such other harmonics may tend show similar behaviour as IER changes. In some embodiments, discerned characteristics 120 other than magnitude may be used to determine the ratio (e.g. the integral over an interval near a peak or over a neighbourhood, as described above).

[0105] Of course, the definition of the ratio as Pi Po is just one possible formulation. The ratio may be defined as P ₀ /Pi (i.e. inverted) or as some other relationship between characteristics of the frequency domain data 116 in the fundamental and in the harmonics, which in some embodiments may be characteristics other than magnitude. Without loss of generality, and for the convenience of the following disclosure, it is assumed the ratio Pi/Po is being used. Where a different ratio is used, the following disclosure can be varied as appropriate. For example, if the ratio is inverted (e.g. if the ratio is defined as P0/P1), then the relationship between the ratio and ATER estimate 124 is reversed so that a positive change in the ratio corresponds to a positive change in the IER, and a negative change in the ratio corresponds to a negative change in the IER.

[0106] Accordingly, block 122 may determine ATER estimate 124 based on a change in a ratio of some or all of the discerned characteristics 120. For example, block 122 may determine (or receive as part of discerned characteristics 120) the ratio of the magnitude of a first set of frequency domain data 116 at the first harmonic frequency to the magnitude of the first set of frequency domain data 116 at the fundamental frequency. That is, block 122 may determine the ratio = Pi/Po for a given period of time represented by the first set of frequency domain data 116. Block 122 may also determine the ratio r ₂ = Pi/Po for a different period of time represented by a second set of frequency domain data 116. For example, block 122 may receive a first set of discerned characteristics 120 based on the first set of frequency domain data 116 and a second set of discerned characteristics 120 based on the second set of frequency domain data 116. Block 122 may determine ratio ri based on the first set of discerned characteristics 120 and ratio r ₂ based on the second set of discerned characteristics 120.

[0107] Without loss of generality, assume that r ₂ relates to a time t ₂ that is later than the time ti to which ri relates. The change in IER between times ti and t ₂ may be estimated based on a change between r ₂ and ri. For example, AIER may be estimated based on a difference between r ₂ and ri (e.g. AYER = g(r ₂ — r- for some function g). Although the foregoing refers to "magnitude" and the corresponding symbols Po and Pi, other characteristics may be used in determining the ratios and r ₂ . For example, the power of frequency domain data 116 may be such a characteristic in embodiments where power is not equivalent to magnitude (e.g. where power is determined by integrating, summing, or averaging over multiple frequencies, such as over an interval around a peak or over a neighbourhood). Determining the power of frequency domain data 116 is discussed in greater detail below.

[0108] For example, a change of—0.05 between ratios ri and r ₂ (i.e. r ₂ — v =—0.05) can be understood, in the exemplary data of Figure 4A, to correspond to a change of +0.2 in the IER of the subject (i.e. AIER = 0.2). Other relationships between changes in the fundamental characteristic(s) and the harmonic characteristic(s) and changes in IER may be used, depending (for example) on the representation of frequency domain data 116, further experimental data, and/or variation between subjects.

[0109] In some embodiments, the ratio may be determined based on a normalized characteristic of frequency domain data 116 at a harmonic frequency, e.g. as shown in Figure 4A. For example, the ratio may be determined based on the normalized magnitude of frequency domain data 116 at the first harmonic frequency (e.g. normalized relative to the magnitude of frequency domain data 116 at the fundamental frequency).

[0110] In some embodiments, block 122 may estimate AIER based on a ratio average. The ratio average may be based on an average of multiple ratios taken from multiple sets of frequency domain data 116. The ratio average may be determined based on previously- determined ratios, ratio averages, and/or discerned characteristics 120. For example, the ratio average may be based on a moving average of previously-determined ratios, the mean of all previously-determined ratios, the mean of a selection of previously-determined ratios, and/or a weighted average of one or more previously-determined ratio averages and subsequently- determined ratios.

[0111] Persons skilled in the art will appreciate that, depending on how the ratio average, ratios, and discerned characteristics 120 are defined with respect to each other in a particular embodiment, the particular embodiment may be determine the ratio average with respect to any of ratio averages, ratios, and discerned characteristics 120, and/or a combination thereof. Accordingly, although the foregoing and following disclosure refers to ratio averages being determined based on previously-generated ratios, it should be understood that other values (e.g. discerned characteristics 120) may be substituted without affecting the result. For example, block 122 may find an average of discerned characteristics and determine a ratio average based on that, without necessarily determining the ratios of the discerned characteristics 120 as an intermediate step.

[0112] Figure 6 shows an example data set in which the ratio average is a moving average. A ratio rj may be determined based on the discerned characteristics 120 generated from temporal window Wj. More particularly, in this example ratio η is determined based on the discerned characteristics 120 generated from the frequency domain data 116 which is a transformation of time domain measurement data 114 associated with temporal window W,. In the example of Figure 6, the ratio average is a moving average comprising the average of the current window W _n and two previous windows W _n-1 and W _n-2 . For example, ratio average ao is an average (e.g. the arithmetic mean) of ratios r ₀ , ri and r ₂ , which are based on the frequency domain data 116 of windows W ₀ , Wi and W ₂ , respectively. Similarly, ratio average ai is the average of ratios r _ls r ₂ , and ¾ and so on.

[0113] In some embodiments, the moving average may be based on a different number of windows. For example, the moving average may be based on a five-window sequence. For instance, to borrow from the example of Figure 6, ao could be determined after window W ₄ is obtained and could be determined based on r ₀ , , r ₂ , r ₃ , and r ₄ . In some embodiments, the ratio average may be based on a non-sequential set of windows. For example, the ratio average may be based on alternating windows (e.g. Wo, W ₂ , W ₄ , etc). In some embodiments, the ratio average may be based on an average other than (or in addition to) the arithmetic mean. For example, the ratio average may be based on a weighted average. For instance, the ratio average may be weighted in favour of more recent windows, so that window W _n has a larger effect on the ratio average than, e.g., window W _n-2 .

[0114] In some embodiments, the ratio average is used in place of (or in addition to) the ratios in the determination of AIER. For example, the change in IER may be determined to be positive if the change in the ratio average is negative. For instance, turning again to Figure 6, if a ₂ <ai then the change in IER between W3 and W ₄ (or, alternatively, between the time period represented by W ₁ -W3 and the time period represented by W ₂ -W ₄ ) may be determined to be positive. More generally, as described above in relation to ratios r ₀ , ri, ..., if a; is less than than a _j (for some i>j) then then IER may be determined to be increasing (i.e. AIER may be estimated to be positive) and vice-versa. This correspondence assumes that the definition of the ratios provided above is being used (i.e. P ₁ P ₀ or a similar definition which grows as the peak at the fundamental shrinks). As noted above, different formations of the ratio and/or ratio average may change the relationship between the ratio average and ATER.

[0115] In some embodiments, block 122 omits from the estimation of AIER changes in the ratio and/or ratio average which exceed a threshold. The basis for this omission is that movement of the subject can, in some circumstances, induce relatively large changes in the ratio and/or ratio average. If the change in the ratio and/or ratio average is so large that it is unlikely or infeasible for the change to be caused by a change in IER, then that change in the ratio and/or ratio average may be omitted. Without loss of generality, this omission process will be described with reference to the ratio average, although it may also (or alternatively) be applied to ratios and/or characteristics, unless the context provides otherwise.

[0116] In some embodiments, block 122 provides a change threshold d ₀ . If the change in the ratio average (or, for brevity, the "change") exceeds do, then one or more ratio averages may be omitted from the estimation of AIER (e.g. because such changes may be interpreted as relating to subject motion rather than a change in IER). In some embodiments, the ratio average corresponding to the later time is omitted. In some embodiments, the changes between ratio averages may be omitted from the estimation of AIER (and not necessarily the ratio averages themselves). For example, if |aj— a _j [ > d ₀ , for some i > j, then aj and/or the change [ a _I — a _j j may be omitted from the estimation of AIER for the time period

corresponding to index i (e.g. window W;). AIER may, for example, be estimated based on the most recent non-omitted ratio averages. For instance, from the foregoing example, AIER may be estimated based on the ratio averages ¾_ι and a _j (or, if a^— 1 = a,-, then a; _→ may be used, for example), provided that those ratio averages have not been omitted. In some embodiments, a new ATKR estimate 124 is not generated when a ratio average is omitted; the previously-generated AIER estimate 124 may optionally be provided as output by block 122.

[0117] Inflammation giving rise to a change in IER tends to develop relatively slowly - e.g. over the course of several minutes or hours. Accordingly, in some embodiments the averaging and/or smoothing processes described herein (including, for example, the determination of ratio averages) may take place over a similarly long period of time. For example, ratio averages may be determined over a period of 30 minutes or more. In some embodiments, ratio averages may overlap - that is, two ratio averages ai and a ₂ may each include a ratio ri in their averages. In some embodiments, ratio averages are non-overlapping. For example, ratio averages may be determined over consecutive periods of time. For instance, ratio average & _\ may be the average of ratios determined over a 30 minute period, and ratio average a ₂ may be the average of ratios determined over the following 30 minute period.

[0118] In some embodiments, block 122 identifies one or more ratios which gave rise to the change exceeding the change threshold d ₀ . Block 122 may determine one or more substitute ratio averages which omit the one or more identified ratios. In some embodiments, the substitute ratio averages replace corresponding omitted ratio averages and may be used in place of those omitted ratio averages for subsequent tests against change threshold d ₀ and/or for the block 122 estimation of AIER. In some embodiments, the substitute ratio averages may be used in place of the one or more ratio averages for the block 122 estimation of AIER, and the one or more ratio averages may be retained and used for subsequent tests against change threshold do.

[0119] In embodiments where omitted average ratios are retained for subsequent tests against change threshold d ₀ , the moving average may be able to "catch up" (over several iterations) to changes with long-term effects on the ratio average. For example, if the subject moves, the ratio between the fundamental magnitude and the first harmonic magnitude may change significantly; when the subject stops moving, all or part of this change may persist.

[0120] For example, consider the following example data set where the ratios are as follows: r ₀ = 0.40, n = 0.42, r ₂ = 0.62, r ₃ = 0.60, r ₄ = 0.64.

This example data illustrates a scenario where a movement of the subject during temporal window W ₂ corresponds to a sharp increase in ratio r ₂ relative to the preceding ratio ri. The sequence of ratio averages ao = 0.41, ai = 0.52, a ₂ = 0.61, a3 = 0.62 may be determined from the example sequence of ratios r _! -r ₄ ; in this example, each ratio average ¾ is the arithmetic mean of two consecutive ratios rj.i and r;.

[0121] Continuing the previous example, if we assume that do = 0.05, then the changes between ao and ai (0.11) and between ai and a ₂ (0.9) may be omitted. Accordingly, AIER estimate 124 may remain unchanged from the value it had prior to the determination of ao (if any). The change between a ₂ and a ₃ , however, is 0.01, which is less than change threshold do. Accordingly, the change between a ₂ and a ₃ may be used to determine AIER estimate 124. Thus, in this example, method 110 may estimate a change in IER corresponding to a change in the ratio average of 0.1 {i.e. based on the change determined during window W ₄ ), but may not estimate a change in IER during windows W ₂ or W ₃ (or may provide the same AIER estimate 124 as had been previously estimated, if any). In effect, in this example method 110 omits the large short-term change caused by the subject's movement (which resulted in a change in ratio averages exceeding do) before "catching up" to a new longer-term baseline following the subject's movement.

[0122] As described above, in some embodiments, method 110 does not provide an AIER estimate 124 while change threshold do is exceeded. Alternatively, or in addition, method 110 may provide the same AIER estimate 124 as was previously provided (if any) while change threshold d ₀ is exceeded.

[0123] As a related example, consider the following example data set with the same ratios as were provided above:

Window (WO Wo Wi W ₂ w ₃ W ₄

Ratio (r 0.4 0.42 0.62 0.60 0.64

Ratio Average (a 0.41 0.52 0.61 0.62

Change Relative to Previous Ratio Average 0.11 0.9 0.1 (|a _f - ai_ _! 1)

This example data is similar to that shown above for the previous example, but it also includes substitute ratio averages. As noted above, method 110 may determine a substitute ratio average ai which omits r ₂ because the change between ao and ai exceeds do (as described above). Substitute ratio average ai may be used in place of ai when estimating the change in IER. As described above, ai may be the average of the remaining ratios (in this case, ri), the average of ri and r _o , or some other value. In this example, assume that ai is the average of the remaining ratios, so that ai = r = 0.42. Thus, the change in ratio averages between a ₀ and ai is 0.01. Accordingly, method 110 may determine A BR estimate 124 based on a change in ratio averages of 0.01 during window W ₂ .

[0124] As noted above, either of ai and ai may be used for subsequent tests against the change threshold do. Assuming that the original ratio average ai is used for such subsequent tests, the next test would be |a ₂ — a _x | = 0.9 < 0.05 = d ₀ . Since this change exceeds change threshold do, method 110 may then determine a different substitute ratio average a". In this example, a" omits ri (since, in this comparison, ri is significantly smaller than ratios r ₂ or ¾). Assuming, without loss of generality, that a" = r ₂ = 0.62, the change in ratio averages between a" and a ₂ is -0.01. Accordingly, method 110 may determine AIER estimate 124 based on a change in ratio averages of -0.01 during window W ₃ . The AIER estimate 124 determined during temporal window W ₄ is the same as the above example, since no further changes exceed change threshold d ₀ and so no substitute ratio averages are required.

[0125] In some embodiments, block 122 determines AIER estimate 124 based on a linear regression of the ratios and/or ratio averages. For example, block 122 may determine a linear regression of all ratios determined over 30 minutes and use the slope of the linear regression to estimate AIER. In determining a linear regression, block 122 may, for example, omit ratios and/or ratio averages as described above. Block 122 may also, or alternatively, include substitute ratio averages, as described above. In some embodiments, block 122 determines a ratio average by determining a linear regression of ratios.

[0126] In some embodiments, block 122 provides a warning threshold di. If the change in the ratio average exceeds di but does not exceed do then sensing system 14 may provide a warning. For example, if d ₀ < |aj— a _j |≤ d for some i and j, then a warning may be provided; this may be referred to herein as the change in the ratio average being in the warning range. Such a change in the ratio average may be interpreted as being a significant change in IER, in which case it may be desirable to alert a user. Block 122 may otherwise behave as described generally herein; in particular, AIER estimate 124 may be determined based on the change |a _; — a _j |. In some embodiments, a warning is only provided if AIER estimate 124 is in the warning range for a certain duration. For example, block 122 may provide a temporal threshold d _t and may provide a warning if the AIER estimate 124 is in the warning range for at least d _t temporal windows. In some embodiments, a warning is only provided if AIER estimate 124 is negative. In some embodiments, sensing system 14 and/or another sensing system may provide a warning based on a detected duration of the subject's inspiration/inhalation in combination with (or instead of) AIER estimate 124. For example, an unusually long inspiration time coupled with a decrease in IER may suggest a restriction of the airways of an asthmatic patient, and so block 122 may provide an appropriate warning if such circumstances are detected and/or estimated.

[0127] In some embodiments, ACER is estimated based on opposing changes in the fundamental and one or more harmonics. For example, AIER may be estimated based on a change in the difference between the magnitudes at the fundamental frequency and one or more harmonic frequencies of frequency domain data 116. That is, the difference P ₀ — P _x may change between times ti and t ₂ . An increase in that difference (i.e. a positive change) may correspond to an increase in IER (i.e. a positive value of AIER). A decrease in that difference (i.e. a negative change) may correspond to a decrease in IER (i.e. a negative value of AIER). Persons skilled in the art will understand that the converse statements apply to the difference P — P ₀ . Because the characteristics of frequency domain data 116 may be sensitive to movement, interference or other noise, these difference-based estimations may not be as robust a measure of AIER as the ratio-based estimation described above.

[0128] As discussed above, the discerned characteristics 120 may comprise the magnitudes of the fundamental and/or harmonic(s) in frequency domain data 116 (e.g. in a PSD estimate). In some embodiments, block 122 determines AIER estimate 124 based on changes in the magnitudes of the fundamental and the first harmonic (i.e. P ₀ and Pi) in a PSD estimate. Depending on how sensing system 14 identifies peaks, P _o and Pi may be the greatest magnitudes of their associated peaks, and/or they may be the magnitudes corresponding to the particular frequencies identified as the fundamental and/or harmonic frequencies (i.e. R o and RRi).

[0129] Other discerned characteristics 120 may be used instead of, or in addition to, magnitude in the systems, methods and examples described above. In some embodiments, the discerned characteristics 120 comprise powers of the peaks associated with the fundamental and/or the harmonics. In some embodiments, as mentioned above, the discerned

characteristics 120 comprise powers of neighbourhoods of the peaks.

[0130] As will be understood by those skilled in the art, where frequency domain data 116 comprises a PSD estimate (and/or where frequency domain data 116 may be analyzed in order to determine an estimate of the power, as is the case with the DFT), the magnitude at a particular frequency may be understood to correspond to the power of that frequency.

However, as shown in Figure 3, a peak in the frequency domain data 116 may be spread across a frequency range. Accordingly, in some embodiments, block 122 may determine a power for each peak of interest (e.g. the peaks surrounding RRo and RRi). Such an estimate may, for example, be based on an estimate of the integral of the PSD under a peak in the vicinity of RR ₀ and RRi, as described in greater detail below. In some embodiments, block 122 may determine a power of the frequency domain data 116 based on neighbourhoods of the fundamental and/or harmonic frequencies (without necessarily detecting peaks and/or specifically identifying the power of a peak), as described in greater detail below.

[0131] In some embodiments, the power of a peak may be estimated by, for example, summing some or all of the magnitudes of PSD samples in a (frequency domain) vicinity of the peak, and/or some or all of the magnitudes of PSD samples in a (frequency domain) vicinity of the peak exceeding a power threshold. Returning briefly to the exemplary Figure 3, the power of the peak around frequency RRo (i.e. the fundamental) may be estimated, for example, by summing each of the frequency domain samples between sample 170 and sample 172 (inclusive), which exceed exemplary power threshold 140.

[0132] A suitable power threshold may be determined as described above, with reference to identification thresholds 132, 134, 136 (e.g. based on an average, relative to a particular magnitude, user-configurable, etc.). A power threshold may be the same as or different than identification thresholds 132, 134, 136; for example, in some embodiments identification threshold 136 may be used to identify peaks, and a power threshold equivalent to identification threshold 134 may be used to determine the power of a peak. In some embodiments, a power threshold (e.g. exemplary power threshold 146) may be set at or above a noise floor of the frequency domain data 116 so that noise in frequency domain data 116 does not excessively impact the power estimation.

[0133] In some embodiments, different power thresholds may be used for identifying the powers of different peaks. For example, the power threshold used to identify the power of a peak may be determined based on a magnitude P of the peak (e.g. based on the largest magnitude of the PSD samples in or around the peak). In some embodiments, the power threshold used for identifying the power of a peak may proportionate to a magnitude of the peak. For example, in Figure 3, exemplary power threshold 140 for determining the power of the peak surrounding the fundamental may be set by determining ½P ₀ (i.e. one-half of the largest magnitude of the peak). Exemplary power threshold 144 for determining the power of the peak surrounding the first harmonic may similarly be set by determining ½Pi. Of course, other proportions (such as 20%, 25%, 30%, 1/3, 60%, 2/3, 75%, and/or any other suitable proportion) may be used.

[0134] As another example, a power threshold may be determined based on an offset from a magnitude of the peak. For instance, exemplary power threshold 142 for determining the power of the peak surrounding the fundamental may be set by determining Ρο-Δ for some offset Δ; in the example shown in Figure 3, Δ is 6xl0 ⁹ in the units of the vertical axis. Δ may be any suitable value; for instance, in a particular example embodiment in which frequency domain data 116 comprises a PSD, Δ may (for example) be 6 dB. Different offsets may be used for different peaks; for example, an offset used to establish exemplary power threshold 144 for the peak surrounding the first harmonic may be less than an offset used to establish exemplary power threshold 142. For instance, returning to the previous example, the offset for determining power threshold 144 may be 2 dB. Of course, offsets may take values other than 2 dB and/or 6 dB.

[0135] In some embodiments, the power of a peak may be estimated by estimating the integral of the frequency domain data 116 around the peak. This may, for example, comprise summing the magnitudes of PSD samples in a vicinity of the peak, as described above.

Alternatively, or in addition, estimating the integral may comprise estimating the integral of the frequency domain data 116 based on an interpolation of frequency domain data 116. For example, block 122 may estimate the integral of the interpolation of frequency domain data 116 shown in Figure 3 by estimating the integral of the interpolation between the frequencies associated with points 180 and 182 (which, in this example, are the points where the interpolated data intersects with identification threshold 132). Such interpolation and integration may, for example, be performed or estimated according to any applicable methods of numerical analysis.

[0136] In some embodiments, the power of a peak may be estimated by determining an average of one or more of the magnitudes of the PSD samples in a vicinity of the peak. For example, the power of a peak may be estimated by determining the average of the magnitudes of the PSD sample with the lowest frequency and with a magnitude greater than the threshold (e.g. sample 170 in Figure 3) and the PSD sample with the highest frequency and with a magnitude greater than the threshold (e.g. sample 172 in Figure 3) together with the PSD sample having the largest magnitude between the lowest and highest frequencies (e.g. Po in Figure 3).

[0137] In some embodiments, the power of a neighbourhood of a frequency may be estimated as described above, but over a larger or smaller neighbourhood in the frequency domain. For example, turned to Figure 3, block 122 may identify a boundary frequency 190 that splits the domain of frequency domain data 116 into two neighbourhoods - a

fundamental neighbourhood (comprising all frequencies less than boundary frequency 190) and a harmonic neighbourhood (comprising all frequencies greater than boundary frequency 190). Boundary frequency 190 may be included in both, one, or neither of the fundamental and harmonic neighbourhoods. The example boundary frequency 190 was determined by identifying the fundamental frequency / and setting boundary frequency 190 to be 1.5/ (thus placing it roughly equidistant between the fundamental and the first harmonic).

[0138] Other boundary frequencies 190 may be used - e.g. boundary frequency 190 may be pre-determined, determined based on the frequencies of the fundamental and/or harmonics, and/or in other ways. In some embodiments, more than two neighbourhoods are used (e.g. each harmonic may have its own associated neighbourhood). In some embodiments, the union of the neighbourhoods encompasses the entire domain of frequency domain data 116. In some embodiments, the union of the neighbourhoods does not cover at least some of the domain of frequency domain data 116. Neighbourhoods may be overlapping or non- overlapping, depending on the embodiment. [0139] As described above, the power of the signal may shift from the fundamental to harmonic components of the signal as IER changes. Accordingly, discerned characteristics 120 may comprise the power of one or more fundamental neighbourhood(s) and the power of one or more harmonic neighbourhood(s), and these discerned characteristics 120 may be used to determine a change in IER substantially as described herein. The power of a

neighbourhood may be estimated in any of the ways described herein (e.g. via integration, numerical analysis, averaging, summing, and/or otherwise).

[0140] As described above, one or more power thresholds (e.g. exemplary power threshold 146) may be applied to frequency domain data 116. In some embodiments, only magnitudes of frequency domain data 116 which exceed the power threshold(s) are used to determine the power of a neighbourhood.

[0141] In some embodiments, spectra may be averaged over time (e.g. over several spectra associated with several corresponding time domain windows) to smooth the spectra and the resultant averages could be used to obtain the power of peaks according to any of the other techniques described herein.

[0142] Although the preceding disclosure refers to the PSD and magnitudes thereof, the power of a peak may be estimated where frequency domain data 116 is not a PSD estimate. Although the term "power" is used for convenience, power estimated by block 122 is not necessarily analogous to or representative of power in a physical sense (e.g. signal power), and should be broadly understood to refer to a magnitude or a related quantity (e.g. a multiple of the magnitude, the magnitude raised to a power, a translation or transformation of the magnitude, and/or any other function of the magnitude). For example, the power of a peak in a frequency domain data 116 representing an unmodified DFT may be estimated as described above (e.g. by summing the magnitudes of the peak, integrating, finding an average, etc.), notwithstanding the fact that magnitudes of a DFT are not necessarily a direct representation of signal power.

[0143] Certain implementations of the invention comprise computer processors which execute software instructions which cause the processors to perform a method of the invention. For example, one or more processors in a sensing system or a monitoring system may implement data processing steps in the methods described herein by executing software instructions retrieved from a program memory accessible to the processors. The invention may also be provided in the form of a program product. The program product may comprise any medium which carries a set of computer-readable signals comprising 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, physical (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, or the like. The instructions may be present on the program product in encrypted and/or compressed formats.

[0144] Where a component (e.g. a software module, controller, processor, assembly, device, component, 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.

[0145] While a number of exemplary aspects and embodiments are discussed herein, those of skill in the art will recognize certain modifications, permutations, additions and subcombinations thereof. For example:

• The description set out above refers at various times to particular example processes by which aspects of the respiration rate (e.g. fundamental and harmonic frequencies and their characteristics) of the subject may be determined. The present disclosure is not limited to specific processes for estimating the respiration rate. Any process for determining the respiration rate which may provide discerned characteristics 120 from which AIER estimate 124 may be generated may be used.

• The description set out above refers at various times to magnitude of the frequency domain data 116 (e.g. magnitudes of PSD samples). As noted above, magnitude is just one example of a characteristic of the frequency domain data 116 which may be discerned and used in block 122 to generate AIER estimate 124 - other such characteristics include, for example, power (as defined broadly above) and other characteristics of the signal which may be of interest to those skilled in the art. • In some circumstances, the order of operations described in the examples above may be varied. For example, the magnitude-squared DFT may be determined by applying the DFT to time domain measurement data 114 (and/or windows thereof), as described above, squaring the magnitude of the result, and output the magnitude- squared DFT as frequency domain data 116 (perhaps with some additional processing, as described above). However, in some embodiments, sensing system 14 may apply the DFT to time domain measurement data 114 (and, optionally, may take the absolute value of the result) to obtain frequency domain data 116. Such embodiments may still, optionally, determine the change in IER based on characteristics of the magnitude-squared DFT; for example, such embodiments may identify magnitudes of interest in frequency domain data 116 and, during the application of method 110 and/or method 200, may determine the square of those magnitudes without necessarily modifying frequency domain data 116 to represent the magnitude-squared DFT. In effect, in some embodiments, the magnitude may be squared on an ad-hoc basis without necessarily requiring that frequency domain data 116 be represented in that form.

[0146] While a number of exemplary aspects and embodiments have been discussed above, those of skill in the art will recognize certain modifications, permutations, additions and subcombinations thereof. It is therefore intended that the following appended claims and claims hereafter introduced are interpreted to include all such modifications, permutations, additions and sub-combinations as are within their true spirit and scope.