SENSING BLOOD FLOW AND HEMODYNAMIC PARAMETERS

[0001] This relates to apparatus and methods for measuring, monitoring and analyzing biological and physiological functions of the human body, and to biomedical

instrumentation and signal processing relating thereto.

BACKGROUND

[0002] The usefulness of hemodynamic data challenges the art to find methods and devices for obtaining, isolating, determining and monitoring such data more simply, economically and more efficiently. Hemodynamics as discussed herein includes the study of blood flow-related data directly or indirectly related to blood flow, such as: heart stroke volume, cardiac output, pre-ejection period, contractility (ability of heart to contract, inotropy), and related causal or caused bodily dynamics such as exercise and exercise recovery, and the Valsalva maneuver (such as when pushing or straining while holding one's breath, or otherwise doing the maneuver in a medical test). Measurement of blood flow, hemodynamics and cardiovascular performance is integral to a full assessment of an individual's health. Specifically, patients with past conditions of heart disease like heart failure (potentially arising out of one or more of many causes like coronary artery disease, heart valve or heart muscle disorders, past myocardial infarction, hypertension etc.) may need constant monitoring in order to improve a person's quality of life via timely and appropriate diagnostic interventions. While the physiological mechanisms underlying these conditions are fairly well understood, the technology to monitor these physiological vitals calls for considerable improvement. Impedance cardiography, doppler echo cardiography, continuous blood pressure monitoring, etc., are believed obtrusive and expensive, and poorly- suited for ambulatory monitoring applications.

SUMMARY

[0003] Generally, and in one form of the invention, an electronic monitoring device includes an electronic processor having at least one signal input for body monitoring. A memory is coupled to the electronic processor and provides instructions for the electronic processor to isolate a cardiac signal including cardiac pulses combined with other cardiac signal variations. The electronic processor is further operable to execute a filter that separates a varying blood flow signal from the cardiac pulses, and to output information based on at least the varying blood flow signal.

[0004] Other devices, sensor assemblies, electronic circuit units, and processes are also disclosed and claimed.

BRIEF DESCRIPTION OF THE DRAWINGS

[0005] FIG. 1 is a partially-block, partially-pictorial, partially graphical depiction of an inventive structure and process for separating a blood flow signal from heart and other signals using one or more axes of a single chest accelerometer sensor; and a time- magnified waveform inset is a voltage- versus-time graph of three concurrent waveforms including an ECG signal, a filtered heart signal from the Z-axis sensor, and a blood flow signal filtered from Y-axis of the accelerometer sensor.

[0006] FIG. 2 is a block diagram of an inventive wireless system structure and process and including inventive structures and processes from the other Figures.

[0007] FIG. 3 is a combined flow diagram of inventive processes for separating a heart signal from body motion and noise using Z-axis accelerometer sensor input and separating a blood flow signal using Y-axis accelerometer sensor input and further electronically processing the heart signal and blood flow signal jointly for

hemodynamics.

[0008] FIG. 4A is a flow diagram of inventive process for separating a heart signal from body motion and noise using Z-axis sensor input such as for use in FIG. 3 or with FIG. 4B.

[0009] FIG. 4B is a flow diagram of inventive process for separating a heart signal and blood flow signal from each other in Y-axis sensor input, for hemodynamics and heart rate.

[0010] FIG. 5 is a block diagram of an inventive structure and process for variably combining accelerometer signals from multiple axes in various proportions. DETAILED DESCRIPTION OF EXAMPLE EMBODIMENTS

[0011] Example embodiments are described in an implementation for estimating blood flow and hemodynamic parameters from a single chest-worn sensor. The described embodiments enable measurement of blood flow trends (stroke volume, cardiac output) and other hemodynamic parameters (contractility, pre-ejection period, iso- volumic contraction interval) in a non-invasive and minimally obtrusive way. In an advantageous approach, a single miniature sensor (such as a MEMS accelerometer coupled with a data acquisition signal processing embodiment) extracts hemodynamic parameters from an in- plane-of-chest accelerometer Y-axis parallel to a superior-inferior line from head-to-feet.

[0012] These benefits are obtained by themselves and with heart pulse detection, respiration detection, motion detection and other features. Ambulatory and other monitoring of cardiac activity including ambulatory hemodynamic monitoring can find widespread applications in home health monitoring of patients with a history of cardiovascular conditions, monitoring older adults, ICU and hospital monitoring, monitoring vital signs in gurneys or rolling patient transfer beds, in mobility aids like scooters and wheelchairs, ambulances, at accident and trauma sites or in a clinic, allied medical testing facility, airport, etc., and for fitness monitoring at exercise centers, stores and elsewhere. Such capability is directly relevant in applications that involve ambulatory monitoring of cardiovascular and cardio-respiratory health. Among advantages of some embodiments, are: a) single sensor and a single signal used to extract several hemodynamic vitals such as any, some or all of changes in stroke volume, changes in cardiac output, heart rate, isovolumic contraction interval, etc.; b) minimal obtrusion (miniature sensor taped on the chest); c) minimal inconvenience in continuous monitoring applications; and d) low cost and convenience with sensors carried by disposable patches and tapes.

[0013] Isovolumic contraction interval (IVCI) is the duration in early systole when the heart ventricles contract with no change in volume. During IVCI, myocardial muscle fibers shorten but have not developed enough pressure in the ventricles to open the aortic and pulmonary valves. IVCI occurs after closure of the mitral and tricuspid valves and before opening of the semi-lunar valves. IVCI is electronically measured or estimated as the time difference between SI peak PI and beginning time Fl of first peak of the vertical Y-axis accelerometer blood flow signal. ICVI is expected to correlate well with the time difference Fl-Pl (see FIG. 1 inset).

[0014] Stroke volume (SV) is the difference between the end diastolic volume and end systolic volume and is a measure of the blood pumped by the heart per cardiac cycle. A conventional pulse contour method calculates a blood flow variable (milliliters/sec) from the pressure signal and computes the stroke volume by integrating the blood flow signal over a cardiac period. By embodiments of structure and process herein, a peak amplitude PAmp and Jamp of the flow signal derived by filtering from the accelerometer sensor is used to compute relative changes in SV. SV is computed by first applying a blood pressure signal, or a signal related thereto, to a model of the arterial system. One such model is a Windkessel model, non-linear to model behavior of the arteries as they expand under blood pressure, which regards the blood pressure as analogous to a voltage applied to a series-parallel network having a series impedance to an output, and a parallel resistor-capacitor combination across the output. Blood flow is analogous to the voltage across the output of the circuit. The integrated output voltage over a period of heart rate S l-to-Sl is related to stroke volume for that period and is repeatedly computed. Some other models analogize a reflective electrical transmission line to the arterial system. Any model appropriate to the purposes at hand is employed.

[0015] Cardiac output (CO) is the volume of blood pumped by the heart per minute, and defined as the product of stroke volume and heart rate (HR).

[0016] Pre-ejection period (PEP) is the time interval between onset of ECG QRS complex and the cardiac ejection. PEP is calculated as R-Fl interval from the beginning of the ECG QRS complex to the beginning of the first peak in the accelerometer signal (see FIG. 1 inset).

[0017] The electronic measurement in one example jointly processes the flow signal from filtering the accelerometer-derived Y-axis (FIG. 1); and the simultaneous ECG as shown in FIG. 2 or estimating from IVCI as described elsewhere herein.

[0018] Ventricular contractility (VC) measures intrinsic ability of the heart to contract and can be estimated from stroke volume (SV). Increase in stroke volume causes an increase in contractility. Contractility (VC) is measured herein by trending the pre- ejection period PEP or estimating from IVCI as described elsewhere herein.

[0019] FIG. 1 illustrates a system embodiment that provides a convenient monitoring device and measurement set-up. A miniature (weight = 0.08 gram, size = 5x5x1.6 mm) triple axis, low-power, analog output MEMS accelerometer (e.g., STMicroelectronics Part No. LIS3L02AL) is taped onto the chest (e.g., a few inches to the left of the sternum along the third or fourth rib). (Taping the accelerometer sensor or using a chest band presses the accelerometer sensor to or against a bare or shaved portion of the chest and efficiently couples chest acceleration to the sensor.) An acceleration signal along the Z- axis— the dorso-ventral direction orthogonal to the plane of the chest— corresponding to the heart sounds is captured. Simultaneously, the acceleration signal along the Y-axis— in the plane of the chest and oriented vertically upwards— corresponding to the blood flow and related hemodynamics is also captured by the same MEMS (micro- electromechanical system) accelerometer. In this way, the acceleration signal along the in-plane vertical axis (Y-axis) is also captured.

[0020] The chest acceleration signals from both axes are, for instance, concurrently AC coupled (high pass roll-off was dropped about lOx in an example compared to a non- critical 3 Hz) and separately and in parallel are amplified with a gain of 100 and low pass filtered - for anti-aliasing - through a three-stage, five-pole Sallen-and-Key Butterworth filters with a 1 kHz corner frequency. Two commercial quad operational amplifier packages (e.g., Linear Technology Part No. LT1014CN) are used for the analog front- end. The accelerometer signals are then each sampled at 10,000 samples/sec using a data acquisition card (e.g., as available from National Instruments, Austin, TX) and captured and stored on a computer using MATLAB software (Version 2007b, The Mathworks, Natick, MA).

[0021] The AC coupling with approximately 3 Hz cutoff, which is a non-critical roll- off frequency, is provided, for example, by a series coupling capacitor C coupled to an input resistance established for the amplifier.

[0022] In an example for detection of primary heart sounds and cardiac activity (illustrated by FIG. 3 (left side) and FIG. 4A), the acceleration signal is digitally low pass filtered in a step 110 at 50 Hz (using, e.g., a 3326 tap digital FIR filter with a steep 80 dB roll-off over 20 Hz) and decimated in a step 120 by a factor of 10. Roll-off frequency less than 60 Hz attenuates 60 cycle USA power line interference with biomedical signals of interest, and roll-off may be made less than 50 Hz for applicable countries using 50 Hz. While the roll-off frequency could be made higher, this finite impulse response (FIR) filter also desirably attenuates white noise above the frequency range of the signals being monitored. Also in a Phase 1, a high order Savitzky-Golay (S- G) polynomial smoothing filter 130, using 28th order and 401 point frame, is used to capture the relatively slow- varying motion wander, and leave out the more rapidly varying heart sound signal components. (MATLAB syntax for such filter is g = sgolayfilt(X,28,401), where g is the filter output and X is a latest input column vector of 401 sample values of windowed data.) In a Phase 2, the smoothing filter 130 output is subtracted in a step 140 from the decimated low pass filter (LPF) output to obtain heart sounds SI and S2. A folded correlation process in a step 160 then enhances and strengthens the polynomial filtered S1/S2 peaks in the motion-removed acceleration signal. Then the location of the peaks is threshold-detected in a step 170 using an electronic amplitude-based peak picking process, and the selected peaks are counted in a step 180 to calculate heart rate (HR).

[0023] A reference electrocardiogram (ECG or EKG) is acquired simultaneously in a three-electrode (single lead) ECG amplifier by analogous signal path for system testing to compare with the accelerometer-derived cardiac signal and independently extract, e.g., heart pulses, pre-ejection period (PEP), contractility (VC).

[0024] The system extracts primary heart sounds (S 1 and S2 produced by the heart valve pairs closing at the ends of the diastolic and systolic periods respectively of the cardiac cycle) using a Z-axis sensor of the accelerometer worn on the chest (see FIG. 3 left side and Appendix). (The term "heart sound" refers in an expansive way to a signal analogous to cardiac SI, S2, and/or heart murmur or other cardiac waveform features, obtained from the processing of accelerometer data or other sensor data, and not necessarily to an audible sound.) Heart sound SI is concurrent with tricuspid and mitral valve activity and shows on a seismocardiogram as a pulse bundle. Heart sound S2 is mainly associated with pulmonary valve and aortic valve activity. [0025] Hemodynamics from Y-axis: Striking differences appear between the signals picked up by the Y-axis (in the chest plane) and the Z-axis. SI and S2 align well in shape and timing in both axes, but the additional features (lower in frequency with the biggest peak between SI and S2) occur very strongly along the Y-axis, but are hardly present in the Z-axis direction, as of a blood flow-related signal. The phrase flow signal or flow- related signal or blood flow signal is a useful identifying label for this chest acceleration component, recognizing that the accelerometer is sensing a component of chest acceleration in meters/sec . This time- varying acceleration component appears to be a body-reaction acceleration or skin-shear acceleration approximately parallel to the superior- inferior body axis in response to blood flow and is related to an overall force of cardiac contraction in newtons and/or to systolic blood pressure in newtons/square meter and/or to blood flow acceleration in milliliters/sec and/or blood flow velocity in milliliters/ sec. Dynamics of the blood varies spatially and with time in the interior and along the blood vessels of the arterial system emanating from the heart. Arterial wall friction and elasticity are involved, and these change if hardening of the arteries occurs. Accelerometer signal(s) are accordingly processed as for FIGS. 3-5 and post-processed and interpreted, e.g., using a system of FIG. 2 to fully realize benefits.

[0026] FIG. 2 shows a wireless system embodiment 600 for a hemodynamic, cardiac, respiration and motion monitoring system with remarkable component embodiments. (Some system embodiments can be wired embodiments.) Accelerometer sensor 610 and its electronic circuit also have a Bluetooth or other short range wireless modem wirelessly coupled to another short range wireless modem 622 coupled via a streaming data and control interface 624 to a data acquisition signal processing unit 620. Further, modems 640.2 and 670 provide wireless transmission and reception so data acquisition signal processing unit 620 communicates via its modem 640.2 to the remote wireless transceiver unit or modem 670 that is coupled to be one or more display units 650.i and their storage unit(s) 660.L In this way, telemedicine applications are supported.

Acquisition signal processing unit 620 and its modem 640.2 are suitably provided in a residence or ambulance or on the person or in a wheelchair or gurney. The wireless transceiver 670 and display unit(s) 650.i are suitably provided in a clinic, hospital, medical monitoring center or otherwise. Either or both ends of the wireless system may be mobile, e.g., a modem 640.3 and alert/processor/display 680 when a professional in a vehicle is urgently needed to review data in case of emergency and to respond with instructions.

[0027] In FIG. 2, devices provide automatic cautionary responses, warnings, and/or automated monitored therapeutic responses. Upon occurrence of undue excursions of one or more measured parameters or relationships among parameters detected by signal processing unit 620, the remote processor 670 alerts any one or more of medical professional, patient, caregiver, and/or family member via a modem 640.3 and alert/processor/display unit 680 by sending an acknowledgable voice call or message and/or written alert or e-mail. Another modem unit 640.1 couples to a telemedicine therapeutic or assistive device 690 for assisting the patient in some pharmacological, informational, or physically assistive way by administering a medicine, or adjusting a dosage or otherwise. In case of excursions that indicate an extreme medical emergency, the data acquisition signal processing unit 620 may be permitted to locally control the therapeutic or assistive device 690 temporarily and in a maximally- safe way until remote commands are received or first responders can arrive. Mere removal or inadvertent accelerometer 610 detachment from chest is distinguished by electronic processing 620 from affirmatively detected excursions. Telecommunications care (tele-care) assistance may be rendered responsive to sensor 620 data by activating motorized apparatus 690 such as to adjust a motorized bed, or move a scooter into proximity for patient use, or servo-mechanically actuate and flush a toilet unit, or open a tub drain to empty a tub, or some other assistance.

[0028] Types of embodiments, among others, can form: 1) a patient- worn medical- sensor and/or therapeutic device that is wired or has a wireless modem; or 2) a patient- worn signal processing and modem module on a belt clip that connects or wirelessly couples to such a sensor and wirelessly couples to a premises gateway or directly transmits to a remote location; or 3) a sensor, signal processor, memory, and modem together in a micro-miniature device that is taped to the chest and communicates to a router or gateway locally, and the latter communicates remotely; or 4) a local router or gateway that includes signal processor, memory, and multiple types of modem to communicate with the sensor and separately communicate remotely, such as in a patient home-based system to telecommunicate to clinic or hospital; or 5) a complete medical clinic system.

[0029] Peripherals and interfaces can include a UART (universal asynchronous receiver/transmitter) data interface and MCSI (multi-channel serial interface) voice and data wireless interface for an off-chip IEEE 802.15 (Bluetooth and low and high rate piconet, Zigbee, and personal network communications) for a wireless circuit. The Bluetooth or Zigbee wireless interface is useful for receiving from and controlling the accelerometer sensor 10 (610) and its associated analog circuitry and digital to analog-to- digital converter ADC in FIG. 1. In arrangements including ECG electrodes and/or a chest microphone, the analog circuitry at the taped-on sensor unit also includes couplings from such pickup elements to the Bluetooth or Zigbee short distance transceiver from the chest sensor.

[0030] Turning to FIG. 3 (right side) and FIG. 4B, the Y-axis signal is filtered or smoothed using a low order (4th order) Savitzky-Golay polynomial filter with a window size roughly 200 milliseconds. Different window length and polynomial orders are feasible. Also, specific polynomials and orders are illustrative and not limiting because they are related to the signal and the sampling rate. The smoothing filter extracts the slow varying (lower frequency) blood flow signal separated from the residue of the heart sounds (SI and S2). Using Savitzky-Golay (S-G) filtering as the smoothing filter is one of various possible ways of extracting the flow signal from the Y-axis. Slow varying respiration and body motion baseline wander (e.g., below about 0.5Hz) is also removed or separated from the blood flow signal in some embodiments by either cascading or combining a high pass filter to attenuate the respiration wander, or using a smoothing filter to isolate the respiration wander and subtracting to yield a residue signal. Further in some embodiments, the respiration is separated from body motion either by detection of variations of SI -SI interval or by monitoring amplitude modulation on the SI peaks.

[0031] An extracted flow signal from filtering of the raw Y-axis signal is called a blood flow component, or flow signal, herein. A 200 msec, window size nicely encompasses either straight or curved portions within a cycle of the oscillating blood flow signal as the raw Y axis signal stream progresses through the filter window. A primary "flow" peak (between S 1 and S2) from the accelerometer-based waveform is prominent on the blood flow signal and bears some relationship to a signal from a ballistocardiogram (BCG) while not being identical to a BCG. This flow signal appears related to acceleration of surge blood as pushed out of left ventricle into aorta. The flow signal has noticeable ringing to it, and can convey some useful information about heart and valve mechanics and possibly other subjects, as described further elsewhere herein. Such information is suitably generated and communicated electronically in the wired or wireless apparatus of FIG. 2 operating according to processes described in FIGS. 3 and FIG. 4A, and/or in FIG. 4B.

[0032] In FIG. 3, a process embodiment is represented physically in system storage unit or memory of FIG. 2 and executed on the signal processing unit of FIG. 2 or digital signal processor DSP of FIG. 1. In FIG. 3 (left side), Z-axis signal processing 110-170 occurs as described in Appendix or elsewhere herein, and additionally the peak detection 170 of the Slpulse is followed in FIG. 3 by a time determination 780 of a first peak location in time PI of the first heart sound peak SI in each heart beat. Counter operation in response to the SI pulses (and S2 pulses or both) electronically generates a heart rate signal as in FIG. 4A.

[0033] In FIG. 3 (right), Y-axis signal processing at step 710 low pass filters (LPFs) the Y-axis input (filter cutoff 50 Hz) and then a step 730 S-G polynomial-filters (or otherwise filters effectively) the LPF signal from step 710 to produce a flow signal output from step 730. A step 770 performs electronic peak (and trough) detection on the flow signal to identify a flow peak amplitude PAmp. A succeeding step 785 identifies the location Fl of that peak in time, as illustrated graphically in FIG. 1 inset. Then, in FIG. 3, both the heart rate signal output HR and first peak location in time PI from step 780 are combined with and/or compared with first flow signal peak time location Fl from step 785 for each heartbeat. The time difference of Fl-Pl represents or is proportional to the isovolumic contraction interval IVCI. IVCI in FIGS. 3 (and FIG. 1) is estimated, for example, as the time difference between time PI of the peak of the SI waveform (derived from the residue signal based on the Z-axis sensor) and the time Fl of the first peak of the flow signal derived from the Y-axis sensor. (In FIG. 4B, some embodiments determine heart sound time PI based on Y-axis processing instead of Z-axis processing.) [0034] In FIG. 3, Stroke volume (SV) and relative changes therein are computed at a step 790 proportional to the peak amplitude PAmp of the flow signal derived by filtering from the Y-axis accelerometer sensor. SV in some embodiments is computed

proportional to the peak amplitude PAmp of the flow signal multiplied by the ICVI estimate. Cardiac output (CO) in FIG. 3 is electronically generated by multiplying stroke volume SV by heart rate (HR) derived from Z-axis processing from the step 140 residue from accelerometer sensor filtering and the result is supplied to display 550. i or 650.i and other devices as discussed in connection with FIG. 2. In FIG. 3, hemodynamic parameters including all of IVCI, SV, CO, PEP, and VC and others are obtained as described herein.

[0035] The inset to FIG. 1 depicts an ECG, heart sound signal, and flow signal from the accelerometer 610. The time PI of occurrence of the SI sound and the peak location in time Fl of the flow signal in FIG. 1 inset are used to provide a time difference for electronically estimating isovolumic contraction interval ICVI. Pre-ejection period PEP is electronically derived using the time difference between the peak of ECG signal and peak Fl of the flow signal. A peak amplitude of the flow signal is signified by PAmp in FIG. 1 inset. Peak amplitude PAmp is electronically derived as the difference between the highest peak of the flow signal during a given heartbeat and a baseline average of the flow signal as indicated by a medial dotted line for that heartbeat. A succeeding peak Jamp is also depicted. An electronic process embodiment generates signals representing a value of PEP for every heartbeat and a value of PAmp for every heartbeat.

[0036] Superimposed flow signal ensembles for individual heart beats aligned to the ECG R-waves show relative changes (up to 2x) in signal amplitude and timing

(jitter/dilation) during recovery of subject from mild exercise. Some embodiments also generate a signal representing the variance of the jitter in either or both of ICVI and PEP during a period of exercise and/or during a period of rest, as well as changes in such variance value, as an indication of heart function and changes therein. In some embodiments that measure pre-ejection period (PEP) and contractility (VC), one of the ECG electrodes is physically combined or associated with the accelerometer sensor of FIG. 1. A miniature microphone can be provided along with the accelerometer in the same chest-worn physical unit 10 for obtaining heart sound audio for parallel processing. Concurrent signals from a single chest accelerometer can be separated to get

physiologically relevant information: heart rate, activity/motion, respiration and intrapleural or intra-thoracic pressure changes, hemodynamics (timings and amplitudes and changes), cough, sneeze, snore, speech, breathing sounds, gait, respiration and other activities.

[0037] Some accelerometer sensor embodiments may be called upon to estimate Pre- ejection period (PEP) and contractility (VC) or trends therein, but lack the ECG electrodes, such at a residence. IVCI is based on S l-to-Fl interval in FIG. 1 inset. IVCI is guardedly either used as a proxy directly or for trends for pre-ejection period (PEP) (and for use in obtaining contractility (VC)) where measurements indicate that it is sufficiently-correlated to PEP or trend therein, or at least sufficient for monitoring variations in PEP and VC to detect excursions from expected parameter ranges needing a clinic visit. Previous clinical measurements with both accelerometer and ECG calibrate a time interval adjustment value a and scale value β to estimate PEP from IVCI according to Equation (15), where measurements indicate that adjustment a and scale value β are statistically significant and sufficiently accurate, and may be measured and downloaded as a table of values a (HR, etc.) and β (HR, etc.) to estimate PEP from IVCI adjusted for heart rate (HR, etc.). Such approach is likely acceptable only when no other alternative is immediately available at the remote location.

PEP = p*IVCI + a (15)

[0038] Some embodiments post-process the peak amplitude PAmp on Z-axis or other- axis signal amplitude (which correlates well with SV and CO) to provide or derive time- varying output signals and displays in FIG. 2. Such signals and displays show estimation for hemodynamic parameters such as SV and CO and others derivable directly at FIG. 3 step 790 from the amplitude/power of the cardiac SI pulse either independently of, or in combination with, information from the blood flow signal. The estimation may differ from SV and CO themselves by an additive constant and a scale factor, for monitoring applications that begin with a pre-existing physiological state and are interested in subsequent variations and/or unusual departures. Notice that the SV and CO

hemodynamic parameters vary much more slowly with time t (e.g., less than 0.2 Hz or less than 0.1 Hz or so) than respiratory variation in the peak amplitude signal PAmp(t), so that SV and CO and respiration can be derived from or filtered out of the peak amplitude signal PAmp(t) and provided for display and recording.

[0039] In FIGS. 4A and 4B, a process embodiment is represented physically in system storage unit 630 or memory of FIG. 2 and executed on the signal processing unit 620 of FIG. 2 or digital signal processor DSP of FIG. 1. In FIG. 4A, Z-axis signal processing 110-180 is generally analogous to Z-axis processing of FIG. 3 (left side), and outputs heart rate at step 180 and need not be further detailed. In FIG. 4B, Y-axis signal processing 910-940, 970, 980 independently also derives heart rate by obtaining a residue signal from the Y-axis input with the S-G filtered signal subtracted out. Electronic peak detection of the residue signal at step 940 is followed by peak detection 970 and counting 980. The heart rate signal output from step 980 is either combined with and/or compared with the heart rate output of FIG. 4A or used instead of and without the heart rate output of FIG. 4A, depending on embodiment.

[0040] Further in FIG. 4B, the S-G filtered signal from step 930 itself is a ringy flow signal interpreted as blood flow and provided as an electronic output 950 for display 650. i and optional storage 660.L Also, flow signal 950 is processed at step 960 to recover a Forcing Function F(t) electronic output indicative of cardiac function. That flow signal is processed at a step 965 to estimate one, two, or all of a triplet of 2 ^nd order model parameters for mass m, dashpot p (rho), spring γ (gamma) as further electronic outputs. FIG. 3 steps 770, 785 and 790 in FIG. 4B use the flow signal from step 930. The outputs can be further post-processed into electronically-represented interpretations and displays in FIG. 2 of internals of chest and heart and states of function. Post-processing is applied to the damped oscillatory flow signal at 950 of FIG. 4B (derived from the Y-axis accelerometer sensor) as from a 2nd order system model that has a forcing function F(t) (newtons) and constant coefficient parameters for mass m, dashpot p (rho), spring γ

(gamma) in its 2nd order linear differential equation of Equation (16). Variable y represents physical displacement of chest sensor from average y position (conceptually measured from some stationary point of reference on the body, such as the hips, relative to which the chest is displaced).

md ² y/d t ² + pdy/dt + y y(t) = F(t) (16) [0041] A model of a standing individual has a first triplet of those parameters subscripted "1"; or lying prone has a second triplet subscripted "2." In standing and prone positions, Y-axis sensor is positioned same on chest, parallel to superior-inferior axis of body symmetry. Flow signal 950, g(i) derived by step 930 from Y-axis sensor (e.g., by S-G filtering), is a series of samples g(t) each substantially proportional to the second derivative d ² x/d t ² itself in Equation (16). Dashpot parameter p introduces energy dissipation, and time constant τ of decay of the damped oscillatory signal is related to mass parameter m (kilograms) divided by the dashpot parameter p

(newtons/(meters/sec)).

T = m/p (17)

[0042] Frequency f _s of the damped signal is related to (1/2π) ₁ | ^" (γ ), i.e., the square root of the ratio of the spring parameter γ (newtons/meter) divided by mass parameter m (kilograms):

f _s = (l/2jc)J ^~ (y/m) (18)

[0043] Post-processing suitably estimates F(t) /m , as by numerical integrations from damped oscillatory flow signal waveform from y-axis, S(t) = d ² y/d t ² using Equation (16) written as Equation (19). The numerical integration begins as each spindle-shaped accelerometer Y-axis waveform commences in FIG. 29 (3rd waveform) for a given heartbeat. Constants of integration are zero position, zero velocity. For applications based on the shape or morphology of the forcing function F(t), the mass m is merely a constant of proportionality that does not affect the shape, or mass can taken as that of head and torso, a fraction (e.g., 0.6) of body mass in kilograms. Time constant τ is numerically estimated as the length of time from peak of the spindle-shaped acceleration waveform to the time when the waveform is about 2/3 dissipated (i.e., reduced on later end of the spindle to 1/e of its earlier peak amplitude, where e is the base of natural logarithms 2.71828...). Time constant τ is about a quarter of a second. The frequency f _s is numerically estimated as the number of cycles in some portion of the spindle-shaped acceleration waveform divided by the time in seconds occupied by that portion. In FIG. 29, the frequency f _s is about 8 Hz.

F(t)/m = S(t) + (l/τ) s{t)dt + (2π f _s f s{t)dt (19) [0044] Forcing function Fy (t) component parallel to the Y-axis sensor may arise by 1) physical acceleration of the heart upon ventricular contraction and 2) the acceleration of blood expelled from the left ventricle surging into the aorta. Parameter γ for spring- constant and parameter p for dashpot relate to some gross average of mechanical properties of the interiors of chest and abdomen. Mass parameter m probably is related or proportional to mass of torso and perhaps head, but probably not to mass of legs because legs not accelerated in the Y-axis direction. SI -SI waveform also has a rising amplitude of oscillation immediately preceding the damped oscillation, see FIG. 1, that can be due to entry of blood from venae cavae into right atrium, and right ventricular contraction into the pulmonary artery. Embodiments for extraction and analysis of forcing function F _Y (t) can provide information on cardiac and pulmonary function and hemodynamic data.

[0045] X-axis sensor information can be monitored as well as Y-axis and Z-axis. An X-axis sensor is oriented to sense acceleration laterally across the chest. A transverse displacement variable x for purposes of Equation (20-X) represents side-to-side physical displacement of the chest sensor from an average x position (or conceptually also from a point of reference such as the center of mass of the heart relative to which the chest is displaced.) Signal from X-axis sensor is filtered in parallel with and analogous to filtering of Y-axis signal. Because of the assymetrical location and slantwise inclination and of the heart in the chest, the filtered X-axis sensor provides a lateral (side-to-side) component F _x (t) of vector forcing function F(t). Forcing function components F _Y (t) and F _x (t) can provide further useful information on cardiac function, pulmonary function, properties of the pleura, pleural cavity, and pericardium, as well as hemodynamics information relating to the aorta, venae cavae, and pulmonary arteries and pulmonary veins by any suitable process. The parameter triplets are respectively subscripted "1Y" and "IX" to designate a standing position ("1") and the Y-axis or X-axis sensor involved. If the prone position is involved then the subscript "1" is changed to "2."

mi _Y d ² y/d t ² + p _lY dy/dt + y _1Y y(t) = Fi _Y (t) ; (20- Y) mix d ² x/d t ² + pix dx/dt + y _lx x(t) = Fix (t). (20-X)

[0046] Turning to FIG. 5, in some embodiments a coin- sized accelerometer sensor 10 and transmitter, transceiver, or transponder chip are light in weight and firmly mounted on a thin, resilient plastic support plate with round smoothed or flanged periphery taped to chest. Accelerometer sensor Z-axis is perpendicular to the plastic support. The chip may harvest power from an interrogation signal, or get power from a small battery also on the plastic support. An MSP430™ processor from Texas Instruments or other low power processor and short distance wireless transmitter can be included and may also have a wireline interface and USB (universal serial bus) connector for instance.

[0047] In FIG. 5, the orientation of the X, Y, and Z axes of the accelerometer sensor on the chest may vary, but merely distributes the overall physical acceleration vector a to the different sensor axes according to their vector components. Accordingly, an electronic processing module 990. i for virtual re- orientation or optimization of the accelerometer sensor signals rotates axes by multiplication by a rotation matrix that combines the signals for X, Y, Z as if the accelerometer axes themselves were rotated, as shown by Equation (21). Angle Θ is an angle by which Y-axis sensor is virtually rotated from its affixed position on the chest to align with the foot-to-head direction on the body. An angle φ is an angle to virtually rotate Z-axis sensor from its affixed position approximately perpendicular to the chest toward that foot-to-head direction. A vector V represents the Z-axis signal, the Y-axis signal and the X-axis signal and is matrix multiplied electronically in FIG. 5 according to the rotation product R*V , using rotation matrix R expressed by Equation (21).

[cos φ (sincp sinG) (sin φ cosG)

_{R =} 0 cos9 -sinG (21)

-sincp (cosc sinB) (coscp cosG)]

[0048] The axis rotations are suitably customized to maximize signal output (e.g., blood flow, heart sounds). The angles Θ and φ are each varied by a respective given feedback control circuit 995. i to maximize the desired type of signal output. In FIG. 5, for instance, some embodiments execute a feedback control 995.1 to thus maximize the heart sound signal for the heart monitoring path, and the axes-rotation parameters for the rotation process 990.1 for a feedback loop 990.1, 130, 140,..., 995.1 established either as a configuration routine before run-time, or dynamically at run-time. The feedback loop rotates the axes to deliver a linear combination of Z-axis and Y-axis (and X-axis can also be useful) as an input in place of the raw Z-axis signal in FIGS. 3 and 4A to the Savitzky- Golay polynomial filter 130 for the Z-axis to maximize the amplitude of the SI peaks at the output of the Folded Correlation, for heart sound and heart rate monitoring purposes. Analogously, in FIG. 5, some embodiments additionally or alternatively execute a feedback loop 990.2, 930,..., 995.2 by independently rotating axes to deliver a linear combination of X-axis and Y-axis (and Z-axis can also be useful) as an input in place of the raw Y-axis signal in FIGS. 3 and 4B to the Savitzky-Golay polynomial filter 930 for the Y-axis, to maximize the blood flow signal peak amplitude PAmp of FIG. 1. One rotation 990 and one feedback control 995 can operate from signals jointly, like heart sounds amplitude and/or blood flow signal amplitude. Modes of operation and configuration can be activated or disabled by one or more control registers with bits or bit fields for them. A manual mode, if activated, overrides feedback controls and lets a clinician manually optimize virtual rotations of signals on the computer display of FIG. 2. FILTERING FOR SIGNAL SMOOTHING

[0049] In FIG. 3 (left side) and FIGS. 4A and 4B, gross motion is tracked and canceled out from the primary accelerometer-based signal to leave behind what is called a residue signal. A polynomial smoothing filter 130 (for example, a Savitzky-Golay filter) electronically and digitally smoothes a given accelerometer-based data signal stream by least- squares approximation within a specified data window by a polynomial of a specified order, and preserves higher order moments around inflection points or at extrema like peaks and troughs that a digital moving average or low-pass filter does not. The residue signal provides primarily the heart sounds and other cardiac activity. Order of about 20th - 30th order is satisfactory; but lower orders and even a lst-order may work for some window sizes and applications. Also, low order polynomial filtering example obtains a blood flow signal for hemodynamic monitoring herein. Using a number of points at least approximately half again (1.5 or more times) an order of the polynomial and even substantially higher than that can help to reduce noise.

[0050] The smoothing filter 130 of FIGS. 3 and 4A or 730 of FIG. 4B is configured for an order M and frame size (number of sample points N). For a fixed window length, Nw in points, a higher order polynomial will fit high frequency components of streaming data better. For a given order M, a shorter window tw of time will allow fitting the high frequency component better. One way of optimization estimates the inherent order of the low-frequency component and picks the smallest window that satisfies Nw > +l and Nw is odd (i.e., Nw = 2N+1). Smaller window size Nw means a smaller number of taps of the multiply- accumulate filter process to do the smoothing filtering. For an accelerometer signal in some applications, order M=l and window size Nw =3 (sampling frequency is 1000 Hz). In some examples herein, higher orders M and window widths Nw are shown.

[0051] In FIGS. 3 (left side) and 4A, motion, heart sounds and heart rate are electronically separated and ascertained from accelerometer 210 data using the following steps: a) low-pass filtering 110 and decimating 120 the accelerometer data; b) Savitzky- Golay filtering 130 to fit the relatively lower frequency motion data; c) subtracting 140 the output of the Savitzky-Golay filter from the low-pass filtered accelerometer data (from step a) to obtain the heart sounds; d) performing 160 folded correlation to enhance the primary heart sounds (SI and S2) peak locations; and e) peak picking 170 to count the number of S 1 peaks in a predetermined or configured segment (time interval) and counting 180 the heart rate HR in beats per minute BPM.

[0052] Decimation passes one sample in every no samples (not necessarily removing all but 1 sample in exactly ten). Thus, if a sample/ ADC delivers f _s samples per second, then a decimation frequency is fs/no samples per second. If a window period is tw seconds, then window points Nw = 2N+1 is Nw = 1+fs tw / no. Select window period tw in view of time period of particular features and behavior of interest to be obtained.

Select sampling frequency fs in view of cost, physical size and complexity of FIG. 1 antialiasing filter AAF at 0.5fs or less, ahead of sampling f _s ). Set sampling frequency f _s greater than the Nyquist frequency for sampling the AAF output. Select decimation parameter no, to yield fs/no higher than 50 Hz LPF after sampling/ ADC for effective LPF. Select decimation parameter no also to yield a number Nw of window points high relative to filter order M to keep filter noise low while Nw window points are still low in number but enough so filter computations still satisfactorily filter the signal stream in the window. Filter computations are related to the product of the number Nw of points per window multiplied by a rate number rw of windows processed per second. If rw = Nw tw, the computations are proportional to Nw / tw , which may motivate fewer window points and longer window times for energy-saving and lower cost. Remarkably, the examples herein satisfy these considerations for some applications and other examples may be devised.

[0053] Processes for various forms of step 130, 730, 930 can be partitioned into offline and real-time online electronic processes. The expression ll(x-Ab)ll in Equation (1) below signifies the sum of squared differences between the [(2N+l)xl] respective data stream vector sample points or stream components and the (2N+1) respective estimates of those stream components provided by multiplying a [(2N+l)xM] transform matrix A times a [Mxl] vector of transform coefficients b _j . The number of transform coefficients b _j is M, and they form a [Mxl] vector b. A gradient v is the [Mxl] vector of first partial derivatives with respect to the transform coefficients b _j . The number M of coefficients b _j is called the order, and if the number of transform coefficients b _j is M, then the order of the process is M. The [MxM] matrix of second partial derivatives with respect to the transform coefficients b _j is signified by v V . The filter procedure involves, and in effect forms, a coefficients change coefficient vector Δ b for updating an initial transform coefficient estimate b=0 (i.e., all coefficients initialized to zero). This procedure pre-multiplies the matrix of second partial derivatives times the negative of the gradient to obtain that transform coefficients change vector:

A b = - (w ll(X-Ab)ll y ^l v ll(X-Ab)ll . (1)

Since the Equation (1) involves a quadratic expression and starts from b=0, the process directly finds the values of the transform coefficients b = Δ b in one pass without iterating additionally.

[0054] Equation (2) represents the result of performing the calculus operations represented by Equation (1). (Some embodiments transmit the coefficients b from Equation (2) to a remote site for record storage and further analysis, since they effectively compress much of the information in the data window. If coefficients are to be transmitted, the [Mx(2N+l)] matrix (A T A) - ^" 1 A T is pre-computed and then multiplied by each data window locally on the fly. Other embodiments omit such compression and/or transmission, or only do it locally on remote command, and thereby save some power and processing complexity.)

b = (A ^T A) ^_1 A ^T X , (2) provided the inverse (A T A) - ^" 1 exists— when the rows of the matrix A are linearly independent (full rank) and enough data points Nw = (2N+1) deliver enough

corresponding columns of the matrix sufficient for an inverse to be delivered.

[0055] In the special case of a polynomial transform process, a matrix of indices is raised to powers, wherein the jth column element A _nj in the nth row of transform matrix A is raised to a power: n ^J . In other words, for the 2N+1 different values of n from -N to +N in the window of a data stream X(i+n), the transform finds a set of coefficients b _j for a well-fitting power series to approximate all the values. Such a power series in general is

X'(i+n) = b ₀ + bi n + b ₂ n ² + b ₃ n ³ + b ₄ n ⁴ +... b _M n ^M . (3)

[0056] Savitzky-Golay filtering outputs as the filter output g(i) for the window indexed i the value of bo estimated by Equation (2) for each data window, and

successively window-by-window for successive indices g(i). Rows of matrix A are orthogonal when the inner product is zero for any pair of different ones of them. These rows are illustrated in Table 1, below. The rows of values A _nj in matrix row n are non- orthogonal for the example of a polynomial transform. (" ^Λ " signifies raising to a power) .

Table 1 : Arrangement of Matrix A

Power m (Order M)

Points 0 1 2 3 ... M

n=-N: [1 (-N) (-Ν) ^Λ 2 (-Ν) ^Λ 3 . , .(-Ν) ^Λ Μ] n=-l: [1 -1 1 -1 . , .. (-1) ^Λ Μ]

n=0: [1 0 0 0 . , , . 0 ^Λ Μ]

n=l: [1 1 1 1 .. . 1 ^Λ Μ ] n=+N: [1 N N ^A 2 N ^A 3 .. . Ν ^Λ Μ].

[0057] Next, the process finds an estimated data stream

X'=Ab = A (A ^T A) ^_1 A ^T X . (4)

[0058] An electronic process is set up in a processing circuit as represented by

Equation (2) and electronically executed by the processing circuit. For Savitzky-Golay filtering, the process is optimized to only find g(i) as the estimated value of bo and also to perform as much off-line pre-computation as possible. Accordingly, Equation (4) is revised as in Equation (5) below, to use only the n=0 row [lxM] of the first pre- multiplying matrix A instead of the whole matrix A in Equation (4),

g(i) = [1 0 0 ...0] (A ^T A)- ¹ A ^T X(i) . (5)

[0059] Sometimes a mathematical presentation of Savitzky-Golay filtering regards the window as multiply- added by a set of (2N+1) filter coefficients c(n). Here, a [lx(2N+l)] filter coefficient vector C is introduced so that g(i) = C X(i) = ∑ c(n)x(i+n), (6)

n=-N

where

C = [1 0 0 ...0] (A ^T A) ^_1 A ^T . (7)

[0060] In Equation (8), an alternative notation CI equivalent to Equation (6) post- multiplies Equation (7) by a [(2N+l)x(2N+l)] identity matrix I and designates each of the columns of that identity matrix I as [(2N+l)xl] unit vectors ε _η . The phrase 'unit vector' for ε _η means a [(2N+l)xl] vector of all zeroes except for a one (1) at the nth row position. Furthermore, only the matrix inversion computations to form the first row of inverse matrix (A T A) - ^" 1 are relevant and are performed, considering the pre-multiplication by [lxM] row n=0 vector [1 0 0...0]. Thus, the filter coefficients are also equivalently expressed in the notation of Equation (8), which is equivalent to Equation (7).

c(n) = {(A ^T A)- ¹ (A ^T ¾ ) }o . (8)

[0061] The Savitzky-Golay filter does a local polynomial fit in a least square sense. For a given input variable data window x(i+n) and window of length 2N+1 and chosen polynomial degree M, the filter output is given by g(i). Filter coefficients c(n)— 2N+1 of them— are computed, e.g. off-line, by electronic operations represented by Equation (7) or (8) and loaded into flash memory of a small signal processing unit either worn on the person or provided nearby and coupled wirelessly to the accelerometer sensor 10. The signal processing unit suitably has a digital signal processor circuit such as processor 220 that electronically performs multiply-accumulates (MACs) represented by Equation (6) according to a stored program accessing the filter coefficients c(n).

[0062] Some other embodiments use windows that are not centered around the value at index n used as the output (e.g., n=0). Other matrix types may be used instead of a polynomial transform matrix A (Table 1). The skilled worker chooses the desired transform, the window (frame) size (e.g., 2N+1) and the order M. Also, note that g(i) output of a first filter procedure produces a data stream that itself can be windowed as represented by column vector gi(i ₂ ) in Equation (9B). Equations (9 A), (9B) cascade two lower order filters of Equation (4) and the electronic processing complexity of the computations in implementation is then minimized. The transform matrices Al and A2 can be the same or different, the window sizes (2Ni+l) and (2N ₂ +1) can be the same or different, and the orders Mi and M ₂ can be the same or different, independently of each other.

gi(ii) = [Ι Ο Ο . , Α χ , αηά (9A) g ₂ (i ₂ ) = [1 0 0 ...0] (A2 ^T A ₂ )- ¹ A ₂ ^T gi(i2). (9B)

[0063] Some embodiments may also apply to the SG process a diagonal weighting matrix W which is all zeroes in a [(2N+l)x(2N+l)] matrix except for weights down the main diagonal. The weights can, for instance, be one at the middle of the diagonal and diminish symmetrically in value farther from the middle of the diagonal. Equation (1) is replaced by:

A b = - (v v llW(X-Ab)ll ¹ v IIW(X-Ab)ll . (10)

[0064] Then the electronic process represented by Equation (5) for the output instead is:

g(i) = [1 0 0 ...0] (A ^T WA) ^-1 A ^T W X(i). (11)

[0065] The selection of transform type and matrix A is fixed/predefined by configuration or determined semi- static manner in some embodiments. Dynamic configuration or selection of the matrix A or transform type or parameters of a given transform is contemplated.

[0066] Some other embodiments store and average a set of values from the transform output of Equation (4) from different windowed segments of the data stream X. Thus, several filters are in parallel and averaged in an offset manner. All the values represent a reconstructed value for a same instant of time (i+n) = t. This approach not only uses the n=0 row to approximate bo but also uses the transform approximations to the other coefficients that are available from Equation (3). The electronic processor 220 (and/or 240) electronically performs Equation (12), where X'(i+n) is from Equation (4). +N

g(t) = [l/(2Ni+l)]∑ X'(i+n)l _{(i+n = t)} . (12)

n=-N

[0067] Thus, other types of processes can accord with the teachings herein, whether Savitzky-Golay or not. The skilled worker sets up a test bench with library

accelerometer-based waveforms and then makes the transform matrix choices, choice of number of points (2N+1), and choice of order- value M, either manually or by an automated process. Filtering choices are tested by visual inspection of a display of output from FIG. 1 process or by automated process according to metrics of false negatives and false positives, etc. Transform matrix A values, and values of (2N+1) and M, 2Ni+l, etc., for one or more such filter processes having favorable metrics are then loaded into the monitoring device flash memory or hard drive and executed in real time on the processor. A transform for an embodiment approximates an actual data stream vector x(i+n) and produces an output signal stream that follows the heart sound peaks well over time in response to a data stream X herein derived from a body- worn accelerometer. Some embodiments have reduced processing complexity by using low enough frame size (2N+1), order M and/or using an efficient transform matrix A to achieve desired performance for the purposes for which the monitoring is intended. The same transform is desirably low-complexity and well-performing over numerous patients, accelerometers and their positioning on the body, and in different environments of use, such as clinic, hospital, home, exercise venue, etc.

[0068] In envelope-based noise rejection 150, an amplitude envelope is generated. The residue signal stream r(i) = x(i)-g(i) of of FIGS. 3 and 4A has some remaining noise from subtracting step 140 that subtracts the smoothing filter output g(i) from the LPF- supplied input x(i). An envelope is fitted to the residue. Amplitude -based processed output R(i) is derived from the envelope-fitted residue r(i). Noise n(i) near the horizontal axis is substantially reduced. Operations can use an envelope-related variable gain function, or generates zero signal output when envelope is below a low threshold that passes peaks.

[0069] In folded correlation 160, an input data stream of residue r(i) = x(i)-g(i) , or envelope-processed residue R(i) comes in, and g(i) is the output of the smoothing filter such as by Equation (6) or from a buffer memory for it. R(i) is windowed by a further data window (a frame) of length 2N ₂ +1 ( with points accessed by an index n=0,l,2... N ₂ ).

[0070] The output f _c (i) of the folded correlation is given by Equation (13):

N2

f _c (i) = ∑ R(i - N ₂ + n) R(i + N ₂ - n) . (13)

n=0

[0071] The digital data stream for heart monitoring residue signal samples R(i) from the smoothing filter subtraction is successively processed in overlapping frames indexed i. The value of 2N ₂ + 1 approximates the width tw of a desired signal event (e.g., an SI, S2, or R-wave). SI is perhaps about 100-150 milliseconds long. If decimated sampling frequency f _s is 1000 Samples/sec, the value of 2N ₂ + 1 is established, e.g., as an odd number between 101 and 151, and N2 is some number between 50 and 75 inclusive. Thus,

N2 = RND(tw fs /2). (14)

[0072] The value of N2 can be configured in flash memory, and can be selected or altered by a local or remote operator using the FIG. 2 system. In the electronic folded correlation process 160 represented by Equation (13), the heart monitoring residue samples R(i + N2 - n) from the latter half of each frame are folded around the center heart monitoring sample R(i) in the frame and multiplied by dot product (sum of products in Eq. (13)) with heart monitoring residue samples R(i - N2 + n) in the earlier half of each frame. The result of the dot product is a folded correlation output signal stream f _c (i) corresponding to instant i of the input residue signal stream R(i) in the center of the frame. Succeeding thresholding passes the SI peaks and counts them. This concludes the Appendix.

[0073] Processing circuitry comprehends digital, analog and mixed signal

(digital/analog) integrated circuits, ASIC circuits, PALs, PLAs, decoders, memories, and programmable and nonprogrammable processors, microcontrollers and other circuitry. Couplings and connections can be ohmic, capacitive, inductive, photonic, and direct or indirect via intervening circuits or otherwise as desirable. Flow diagrams and block diagrams are each interpretable as representing structure and/or process.

[0074] Those skilled in the art will appreciate that many other embodiments and variations are also possible within the scope of the claimed invention. Embodiments having different combinations of one or more of the features or steps described in the context of example embodiments having all or just some of such features or steps are also intended to be covered hereby.