ORIGINAL RESEARCH published: 12 January 2022 doi: 10.3389/fphys.2021.739035 Combining Physiology-Based Modeling and Evolutionary Algorithms for Personalized, Noninvasive Cardiovascular Assessment Based on Electrocardiography and Ballistocardiography Nicholas Mattia Marazzi 1 , Giovanna Guidoboni 1,2*, Mohamed Zaid 1 , Lorenzo Sala 3 , Salman Ahmad 4 , Laurel Despins 5 , Mihail Popescu 6 , Marjorie Skubic 1 and James Keller 1 1 Department of Electrical Engineering and Computer Science, University of Missouri, Columbia, MO, United States, 2 Department of Mathematics, University of Missouri, Columbia, MO, United States, 3 Centre de Recherche Inria Saclay-Ile de France, Palaiseau, France, 4 Department of Surgery, School of Medicine, University of Missouri, Columbia, MO, United States, 5 Sinclair School of Nursing, University of Missouri, Columbia, MO, United States, 6 Department of Health Edited by: Management and Informatics, School of Medicine, University of Missouri, Columbia, MO, United States Kouhyar Tavakolian, University of North Dakota, United States Purpose: This study proposes a novel approach to obtain personalized estimates of Reviewed by: cardiovascular parameters by combining (i) electrocardiography and ballistocardiography Parastoo Dehkordi, Heart Force Medical Inc., Canada for noninvasive cardiovascular monitoring, (ii) a physiology-based mathematical model for Chenxi Yang, predicting personalized cardiovascular variables, and (iii) an evolutionary algorithm (EA) Southeast University, China for searching optimal model parameters. *Correspondence: Giovanna Guidoboni Methods: Electrocardiogram (ECG), ballistocardiogram (BCG), and a total of six
[email protected]blood pressure measurements are recorded on three healthy subjects. The R peaks in the ECG are used to segment the BCG signal into single BCG curves for each Specialty section: This article was submitted to heart beat. The time distance between R peaks is used as an input for a validated Computational Physiology and physiology-based mathematical model that predicts distributions of pressures and Medicine, a section of the journal volumes in the cardiovascular system, along with the associated BCG curve. An EA Frontiers in Physiology is designed to search the generation of parameter values of the cardiovascular model Received: 09 July 2021 that optimizes the match between model-predicted and experimentally-measured BCG Accepted: 14 October 2021 curves. The physiological relevance of the optimal EA solution is evaluated a posteriori Published: 12 January 2022 by comparing the model-predicted blood pressure with a cuff placed on the arm of the Citation: Marazzi NM, Guidoboni G, Zaid M, subjects to measure the blood pressure. Sala L, Ahmad S, Despins L, Popescu M, Skubic M and Keller J Results: The proposed approach successfully captures amplitudes and timings (2022) Combining Physiology-Based of the most prominent peak and valley in the BCG curve, also known as Modeling and Evolutionary Algorithms the J peak and K valley. The values of cardiovascular parameters pertaining to for Personalized, Noninvasive Cardiovascular Assessment Based on ventricular function can be estimated by the EA in a consistent manner when Electrocardiography and the search is performed over five different BCG curves corresponding to five Ballistocardiography. Front. Physiol. 12:739035. different heart-beats of the same subject. Notably, the blood pressure predicted by doi: 10.3389/fphys.2021.739035 the physiology-based model with the personalized parameter values provided by Frontiers in Physiology | www.frontiersin.org 1 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment the EA search exhibits a very good agreement with the cuff-based blood pressure measurement. Conclusion: The combination of EA with physiology-based modeling proved capable of providing personalized estimates of cardiovascular parameters and physiological variables of great interest, such as blood pressure. This novel approach opens the possibility for developing quantitative devices for noninvasive cardiovascular monitoring based on BCG sensing. Keywords: cardiovascular physiology, cardiovascular monitoring, ballistocardiography, physiology-based modeling, evolutionary algorithm (EA), personalized modeling, cuffless blood pressure estimation 1. INTRODUCTION be taken on virtually any platform using electric leads placed on the body in a standard configuration. Recently, a resurgence of Cardiovascular diseases (CVDs) are disorders of the heart and BCG research has occurred, as new sensing devices (e.g., in the blood vessels, including heart failure, stroke and hypertension, form of bed sensors, chair sensors, weighing scales) allow easier, and represent the first leading cause of death worldwide (World noninvasive capture of the BCG signal. Several of these sensors Health Organization, 2021). Early detection and intervention are now available commercially (Alametsä et al., 2008; Chen et al., of CVDs can reduce the number of preventable hospital 2008; Shin et al., 2008; Young et al., 2008; Inan et al., 2009, 2014; readmissions, thereby helping patients maintain a better quality Giovangrandi et al., 2011; Heise et al., 2011; Satu and Jukka, 2012; of life while significantly reducing healthcare costs (Wasfy et al., Zimlichman et al., 2012; Paalasmaa et al., 2014; Helfand et al., 2014; Soucier et al., 2018). 2016; Katz et al., 2016; Huffaker et al., 2018). Cardiovascular function and oxygen delivery to the tissues Deciphering the cardiovascular mechanisms that determine depends on adequate hemoglobin stores, oxygen uptake from the shape of the BCG waveform in a particular individual is the the lungs and cardiac output (CO). This delivery system relies key to fully unlocking the potential of BCG-based monitoring on a complex interplay between the pumping action of the for noninvasive cardiovascular assessment. For example, Etemadi heart and the biomechanical properties of the vasculature et al. showed that the relative time delay between the ECG and (Chang et al., 2002; Vincent and De Backer, 2013). Thus, BCG peaks is an indicator of myocardial contractility (Etemadi effective cardiovascular monitoring should provide a quantitative et al., 2011). Su et al. (2018) showed that the changes in assessment of both cardiac and vascular functions (Holcroft et al., the amplitude of the systolic BCG peaks measured via a 2006; Vincent and De Backer, 2013). Traditional monitoring hydraulic bed sensor correlated with the change in blood pressure techniques, such as electrocardiography, echocardiography and occurring pre- and post-exercise. In this study, we investigate intravascular catheterization, focus primarily on the heart, how to utilize the BCG waveform to estimate cardiovascular providing information on its electrical, mechanical, and fluid- parameters specific to a given subject, which could be used for dynamical functions. A complementary approach is offered an in-depth assessment of cardiovascular function. by ballistocardiography, whose signal, the ballistocardiogram Specifically, we utilize the mathematical model proposed (BCG), captures the repetitive motion of the center of mass of in Guidoboni et al. (2019) to simulate the BCG waveform the human body resulting from the blood motion within the based on fundamental principles of cardiovascular physiology. circulatory system (Starr and Noordergraaf, 1967). Interestingly, The model parameters quantify aspects of cardiovascular the BCG signal reflects the status of the cardiovascular system function that are particularly relevant for CVD monitoring, as a whole, rather than the heart alone, thereby making it such as ventricular elastances and arterial stiffness, and the an ideal complement to traditional monitoring techniques. In model simulations yield predictions of cardiovascular variables, addition, the acquisition of BCG signals is not invasive and such as blood pressure and volumes, and the resulting BCG does not require body contact, thereby eliminating the risk of waveform (Guidoboni, 2020). In Guidoboni et al. (2019), infections and making it a viable option for both hospital and however, model parameters were chosen as representative of in-home monitoring. an idealized individual based on published literature. In this The original device for BCG measurement used by Starr and study, we investigate the use of an evolutionary algorithm (EA) others was a lightweight bed suspended by long cables (Starr to obtain personalized estimates of the cardiovascular model and Noordergraaf, 1967). The blood flow of a subject lying on parameters based on the comparison between model-predicted the suspended bed resulted in the bed swinging; the capture of and experimentally-measured BCG curves on a specific subject. the swing was the BCG signal. A replica of Starr’s suspended Ground truth for the cardiovascular parameters estimated via bed has been built within the MU Center for Eldercare and the EA is not easily available, and this constitutes the major Rehabilitation Technology (CERT) directed by Prof. Skubic. The challenge of this study. Some parameters may be estimated suspended bed is an impractical device for BCG measurement, via noninvasive techniques (e.g., Doppler imaging can be used especially compared to the electrocardiogram (ECG), which can for arterial radii and lengths), which could be utilized on Frontiers in Physiology | www.frontiersin.org 2 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment TABLE 1 | Summary of the participant information involved in this study. Subject Sex Age Weight [kg] Height [cm] 1 Male 25 72.6 189 2 Male 32 72.6 180 3 Male 22 66.2 176 both healthy individuals and patients with CVD. Conversely, some parameters can only be estimated via highly invasive and risky procedures (e.g. ventricular catheterization is needed to assess end-systolic and end-diastolic elastances in the ventricles), which can only be performed when required by the specific health conditions of a patient with CVD. Thus, in order to properly design a study involving multiple techniques to measure cardiovascular parameters on healthy subjects and patients with CVD that could serve as ground truth to validate our EA findings, it is important to assess beforehand which of the many model parameters may be effectively estimated by the FIGURE 1 | The collection Tf M of BCG curves measured experimentally proposed EA method. This constitutes the specific goal of the (gray curves) are reported along with the fkM , with k = 1, . . . , Nc = 5, present investigation. consecutive curves selected randomly as objective curves for the evolutionary The physiological relevance of the personalized solution algorithm (EA) applied to Subject 1. is evaluated a posteriori by comparing the blood pressure estimated via the personalized model with the blood pressure measured via a cuff placed on the arm of the subject. and the high-frequency noise. Cut-off frequencies of 0.7–40 Hz When tested on three healthy individuals, the proposed EA and 1.25–15 Hz have been used for the ECG and BCG signals, performed better in estimating the parameters characterizing respectively (Enayati, 2019). the function of the left ventricle than those of the right The R peaks in the ECG, located via the Pan-Tompkins ventricle. The EA also performed well in estimating arterial algorithm (Pan and Tompkins, 1985), were used to segment stiffness. The satisfactory agreement between model-predicted the BCG signal. Thus, a family of BCG curves is obtained, as and experimentally-measured values of blood pressure support shown in Figure 1. Let us denote by Tf M the family of all the the physiological relevance of these findings and show promise for utilizing the proposed approach as a quantitative method for BCG curves and by f M = f M (t) a single BCG curve in the noninvasive cardiovascular monitoring. family, so that f M ∈ Tf M . The superscript M indicates that these curves are measured, as opposed to those that are predicted by the mathematical model (as shown in section 2.2). We note 2. METHODS that f M has the units of a force [dyne] as it is obtained via the In this section 2, we outline the details on the signal following relationship: acquisition (section 2.1), the physiology-based cardiovascular model (section 2.2), and the design of the EA algorithm f M (t) = m × a(t) [dyne] (1) (section 2.3). where m is the mass of the subject (m = 74.2 Kg for the subject 2.1. Signal Acquisition and Processing for considered in this study) and a(t) is the acceleration [cm/s2 ] ECG, BCG, and Blood Pressure obtained by applying the following conversion to the signal aV (t) Three healthy subjects were recruited for data collection in [V] actually measured by the accelerometer: controlled laboratory settings. The sex, age, weight, and height of the subjects are summarized in Table 1. The subjects were asked aV (t) − offset a(t) = G (2) to lie still on a suspended bed system as previously described S in Guidoboni et al. (2019), while the ECG and BCG signals were recorded. The ECG was acquired via a 3-lead configuration with offset = 2.5 mV, S = 1 mV, and G = 981 cm/s2 (Kionix, and the BCG was acquired with a three-axis accelerometer from 2014). Since the length of a cardiac cycle may vary from beat to Kionix with 1,000 mV/g sensitivity placed on the bed frame of a beat, the length of each BCG curve may not be constant. Thus, suspended bed (Kionix, 2014). The ECG and BCG signals were the k−th curve fkM ∈ Tf M is defined for t ∈ [0, Tck ], where Tck collected simultaneously using an AD Instrument PowerLab is the length of the k−th cardiac cycle computed as the distance 16/35 data acquisition system (ADInstruments, 2014). The ECG between two consecutive R peaks in the ECG. and BCG signals have been filtered via a 6th order Butterworth In addition, the blood pressure was measured via a cuff bandpass filter to remove the low-frequency respiratory motion placed on the arm of the subjects. A total of six blood pressure Frontiers in Physiology | www.frontiersin.org 3 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment measurements were performed, three before the beginning of the arch (i = 3), the thoracic aorta (i = 4), the abdominal aorta ECG and BCG data acquisition and three afterward. An interval (i = 5), the iliac arteries (i = 6), and the cerebral arteries of 5 min was allowed between measurements. It is important (i = 14). The i labels follow those utilized in Guidoboni et al. to emphasize that cuff-based blood pressure measurements may (2019). In this study, for ease of reference, we define the set interfere with the accelerometer-based BCG acquisition, as they I = {2, 3, 4, 5, 6, 14} to indicate all the arteries in the model. introduce spurious movements. Furthermore, it is likely to expect Important parameters for each large artery i, for i ∈ I , are that the blood pressure may decrease as the subject rests on the radius ri and the length li . All arteries are assumed to have the bed for a prolonged period of time. Thus, we adopted a the same Young modulus E characterizing their stiffness. In the protocol for blood pressure measurements to be performed both model, the iliac arteries in the systemic circulation are followed by before and after the ECG-BCG data acquisition on the suspended the resistance R7 representing the peripheral vascular resistance. bed. Ultimately, for each of the six measurements, systolic and Finally, leveraging the electric analogy to fluid flow (Sacco et al., diastolic blood pressures (SBP, DBP) are recorded and the pulse 2019), the model is completed by other capacitors, resistors, pressure (PP) is computed as the difference between the two, and inductors representing the microcirculation and the venous so that PP = SBP−DBP. The average PP value over the six return to the heart. measurements for the same individual is used as a comparison The outputs of the physiology-based closed-loop model with the model prediction, as illustrated in section 3.3. summarized above are the time-dependent distributions of pressures and volumes of blood as it flows in all vascular 2.2. Physiology-Based Cardiovascular compartments. In particular, the computed blood volume in the Model for Prediction of Blood Pressures, left ventricle as a function of time allows us to calculate the Blood Volumes, and BCG Waveform end-diastolic volume (EDV), the end-systolic volume (ESV), the The physiology-based model presented in Guidoboni et al. (2019) stroke volume (SV), the CO, and the ejection fraction (EF). is utilized to simulate the blood flow through the cardiovascular Furthermore, the computed volume waveforms can be used to system and to predict the resulting BCG waveform. In this study, obtain the BCG waveform f P (t) as the acceleration of the center we mention only the features of the model that are relevant for its of mass of the human body resulting from the motion of blood combination with the EA illustrated in section 2.3, directing the through the cardiovascular system. The superscript P indicates interested reader to Guidoboni et al. (2019) for the full details. that this waveform is model-predicted, as opposed to those In the model, the pumping action of the ventricles is described measured experimentally (as shown in section 2.1). Specifically, by the pressure generators f P (t) is computed as UL (t) = ULO aL (t), UR (t) = URO aR (t) (3) X d 2 Vn f P (t) = ρb (t) yn [dyne] (6) dt 2 n∈N where the subscripts L and R indicate the left and right ventricles, ULO and URO are positive constants representing the pressure where ρb is the blood density, the waveforms Vn (t) represent build-up capacity in the ventricles, and aL (t) and aR (t) represent the blood volume occupying the cardiovascular compartment n nondimensional activation functions for the timing of ventricular at time t, with n ∈ N = I ∪ {lv, rv}, and yn represent the contractions defined as distance in the head-to-toe direction between the cardiovascular compartment n, with n ∈ N , and the plane of the heart valves. tanh(qL ta ) − tanh(qL tb ) aL (t) = , We note that N comprises the left and right ventricles, i.e., 2 {lv, rv}, in addition to the large arteries, i.e., I . The values of all tanh(qR ta ) − tanh(qR tb ) aR (t) = (4) model parameters that are not explicitly estimated via the EA 2 described in section 2.3 are assumed to be the same as those for tm = mod(t, Tc ) < Ts and aL (t) = aR (t) = 0, otherwise. In reported in Guidoboni et al. (2019). Equation (4), ta = t − Ta and tb = t − Tb , with Ta , Tb , qL , and qR positive constants. Furthermore, Ts is the length of the systolic 2.3. Evolutionary Algorithm for Searching part of the cardiac cycle, and Tc is the length of the entire cardiac Personalized Parameters in the cycle. Interestingly, Tc can be tailored to a specific individual by Physiology-Based Cardiovascular Model means of the distance between R peaks in the ECG. The pressure An EA is utilized to search for the parameter values of generators defined in Equation (3) are connected in series with the cardiovascular model described in section 2.2 that time-varying elastances defined as yield a satisfactory match between the model-predicted and experimentally-measured BCG waveforms, which we denoted by EL (t) = ELD + ELS aL (t), ER (t) = ERD + ERS aR (t) (5) f P (t) and f M (t), respectively, for a given individual. We recall that an EA is a computational technique that abstracts from where ELD, ELS, ERD, and ERS are positive constants. the mechanisms of evolution to search for optimal solutions to Large arteries in the model play a very important role, as they a problem. The EA search mechanism is inspired by Darwin’s are major contributors to the BCG waveform. Specifically, the Theory of Evolution: similar to individuals evolving inside a model explicitly includes the ascending aorta (i = 2), the aortic species according to their fitness in the environment, solutions Frontiers in Physiology | www.frontiersin.org 4 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment TABLE 2 | Summary of the model parameters in the genotypes considered in the as objective curves for the evolutionary search. An alternative evolutionary algorithm (EA). choice could have been to create a single template representing Genotype Symbol Unit Description the whole curve bundle. We opted for selecting real curves rather than a template in order to preserve the natural curve features Gr ri , i ∈ I cm arterial radii as much as possible, including their different temporal lengths, Gl li , i ∈ I cm arterial lengths while considering multiple curves in order to capture, rather Gy yn , n ∈ N cm coordinates for BCG calculation than discard, beat-to-beat variations. In this study, we considered G lv ELS mmHg cm−3 left ventricular end-systolic elastance Nc = 5. ELD mmHg cm−3 left ventricular end-diastolic elastance It is a physiological fact that no beat is exactly equal to next in ULO mmHg pressure build-up capacity in the left any given person. As a consequence, it is physiological to expect ventricle that part of the genetic code G will present beat-to-beat variations qL s−1 strength of left-ventricular activation and should, therefore, be treated independently when applying Ts s length of the systolic part of the cardiac the EA algorithm to the Nc objective curves selected for the cycle search. This is the case for G lv , G rv , and G vas , which characterize G rv ERS mmHg cm−3 right ventricular end-systolic elastance the strengths of the ventricle contractions and the response of the ERD mmHg cm−3 right ventricular end-diastolic elastance vasculature. The mean radii, lengths, and locations of the major URO mmHg pressure build-up capacity in the right arteries, on the other hand, may vary with age or with the onset ventricle of disease conditions but are not expected to vary over the few qR s−1 strength of right-ventricular activation minutes required for the BCG acquisition. To account for these vas G E mmHg arterial young modulus differences in the expected genotype variations, we proceed as R7 mmHg cm−3 s peripheral resistance follows. We select one of the Nc objective curves and we apply the EA to search for G r , G l , and G y in addition to G lv , G rv , and G vas . The values for g r ∈ G r and g l ∈ G l obtained during this in an EA evolve in the search space in order to optimize an preliminary search are then used to narrow down the interval for objective function (Ranganathan et al., 2018). the radius and length variations considered when applying the The genotype G in our EA comprises the set of parameters EA to the remainder of the curves. In this study, we considered in the cardiovascular model whose optimal values are subject f1M out of fkM , for k = 1, . . . , Nc , for the preliminary search. The to search. In the following, we will denote by g ∈ G a specific information about the height of the subject is used to limit the choice for the parameter values, also referred to as genetic search range for g y ∈ G y in all simulations. The strategy adopted coding. Mutation rules will vary depending on the anatomical in the EA implementation is detailed below. Whenever necessary, and physiological meaning of each parameter in G . To this end, we will highlight the algorithm variations depending on whether it is convenient to write G as the union of the parameter sets it is applied to f1M or to the remainder of the curves. pertaining to: STEP 1: INITIAL POPULATION. An initial population of M = 300 genotypes is obtained by computing mutations according • The radii and lengths of the major arterial segments in the to specific rules. In the following, we will denote by U (a, b) model, defined as G r = {ri , i ∈ I } and G l = {li , i ∈ I }, a uniform distribution in the range [a, b] and by N (0, σ ) a • The coordinates of the cardiovascular compartments with normal distribution with mean equal to zero and SD equal to respect to the heart valves, defined as G y = {yn , n ∈ N }, σ . Uniform distributions will be used in cases where upper and • The function of the left ventricle, defined as G lv = lower bounds for the range of parameter values can be set a {ELS, ELD, ULO, qL , Ts }, priori, either by means of available measurements or domain • The function of the right ventricle, defined as G rv = knowledge. When such information is not available, a normal {ERS, ERD, URO, qR }, distribution is adopted. When the EA is applied to f1M , mutations • The function of the vasculature, defined as G vas = {E, R7 }. are generated according to the following scheme: Finally, we can write G = G r ∪G l ∪G y ∪G lv ∪G rv ∪G vas . A summary of the model parameters in the genotypes that are subject to the gir = U gi,inf r r , min(gi,sup r , gi−1 ) gir ∈ G r i = {3, 4, 5, 6} EA search is provided in Table 2. The phenotype P represents the manifestation of a given (7) genotype which, in our case, is given by the outputs of the cardiovascular model. In the following, we will denote by p ∈ P gil = U gi,inf l l , gi,sup gil ∈ G l i ∈ I , a specific phenotypic instance. The objective function directing (8) the EA search focuses only on the fitness of the BCG waveform y y U (ḡn , α y ḡn ) for α y > 1 ( f P (t), but additional cardiac variables are used to ensure that y gn = y gn ∈ G y , n ∈ N . y y the search results are physiologically acceptable. Specifically, we U (α y ḡn , ḡn ) for α y <1 write P = P fit ∪ P acc , with P fit = {f P (t), t ∈ [0, Tc ]} and (9) P acc = {EDV, ESV, SV, CO, EF}. The overall EA strategy is illustrated in Figure 2. We consider g lv = ḡ lv + N 0, σ lv g lv ∈ G lv Nc consecutive curves fkM (t), with k = 1, . . . , Nc and t ∈ [0, Tck ] (10) Frontiers in Physiology | www.frontiersin.org 5 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment FIGURE 2 | Illustration of how a priori knowledge about the subject, domain knowledge, experimental measurements, and predictions of a physiology-based mathematical model are combined in the EA to identify cardiovascular parameters for a specific subject. g rv = ḡ rv + N 0, σ rv g rv ∈ G rv the arterial tree. Equations (10–12) show that the functional (11) parameters for the ventricles and the vasculature are computed as mutations of the baseline values, indicated with bars, by means g vas = ḡ vas + N 0, σ vas g vas ∈ G vas of normal distributions (as shown in Table 3). In this study, (12) we assumed the SDs to be equal to half of the baseline values. Equation (16) shows that the genetic code g y ∈ G y is selected y via a uniform distribution between the baseline values ḡn , for Equations (7) and (8) show that radii and lengths for each of the y i ∈ N and their scaled valued by means of a factor α , which is major arteries, namely for each i ∈ I , are selected via a uniform assumed to be equal, larger, or less than 1 depending on whether r , g r ] and [g l , g l distribution based on the ranges [gi,inf i,sup i,inf i,sup ] the subject is as tall as, higher, or shorter than a reference height. reported in the literature (refer to Tables 4, 5). An additional In this study, we assumed the reference height to be 180 cm and anatomical constraint has been included in Equation (7) in order we set α y = 1.1 for Subject 1 (height: 189 cm), α y = 1.0 for to ensure that that the radius decreases when moving down Subject 2 (height: 180 cm), and α y = 0.85 for Subject 3 (height: Frontiers in Physiology | www.frontiersin.org 6 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment TABLE 3 | Summary of the baseline values and the SDs utilized in the EA population of M genotypes, denoted as gj ∈ G , with j = 1, . . . , M, simulation. that can be used to start the EA on each of the fkM curves, with Parameter Unit Baseline Value Standard Deviation k = 1, . . . , Nc , individually. STEP 2: PHYSIOLOGICAL CHECK. For each of the genotypes ELS mmHg cm−3 1.375 0.6875 gj in the initial population, the corresponding phenotype pj , ELD mmHg cm−3 0.04 0.02 with j = 1, . . . , M, is computed via the physiology-based ULO mmHg 50 25 cardiovascular model. In order to be physiologically acceptable, qL s−1 6.28 3.14 we require the values of the cardiac variables in P acc to fall within Ts s 0.35 0.175 some broad ranges reported in the literature and summarized ERS mmHg cm−3 0.23 0.115 in Table II of Guidoboni et al. (2019). The rationale behind ERD mmHg cm−3 0.01 0.005 the physiological check is that, in reality, the human body is URO mmHg 26 13 capable of adapting ventricular and vascular parameters so that qR s−1 6.28 3.14 their combined action leads to proper cardiovascular function. E mmHg 3000 1500 Thus, not all randomly selected genotypes may lead to acceptable R7 mmHg cm−3 s 0.35 0.175 results. Even though this step may raise concerns about its applicability in disease conditions, the ranges utilized for the All baseline values are the same as those reported in Guidoboni et al. (2019) except for Ts , whose baseline value has been estimated from the ECG recording. The values for the check are meant to be quite loose and only provide a way SDs are assumed to be half of the baseline values. to exclude obviously erroneous genetic combinations. Thus, as more clinical and experimental data become available on the cardiac variables in P acc in health and disease, they can be TABLE 4 | Summary of the ranges for the arterial radii utilized in this study. used to enrich the domain knowledge and adjust the ranges for Artery Unit pri,inf pri,sup Reference the physiological check. The generation of genotypes for the initial population continues till, after removal of unacceptable Ascending (i = 2) [cm] 1.49 1.91 Wolak et al., 2008 genotypes under the physiological check, a population of M = Aortic arc (i = 3) [cm] 1.14 1.42 Wolak et al., 2008 300 physiologically-acceptable genetic codes is achieved. Thoracic (i = 4) [cm] 0.92 1.28 Joh et al., 2013 STEP 3: FITNESS RANKING. The phenotypes pj are ranked Abdominal (i = 5) [cm] 0.80 1.10 Joh et al., 2013 according to their fitness, which we assume to be the similarity fit Iliac (i = 6) [cm] 0.492 0.725 Joh et al., 2013 between the model-predicted BCG waveform pj ∈ P fit and the Carotid (i = 14) [cm] 0.26 0.36 Krejza et al., 2006 selected objective curve fkM , for k = 1, . . . , Nc . We remind that this is done independently for each of the Nc objective curves. In this study, the similarity is quantified by means of the Euclidean TABLE 5 | Summary of the ranges for the arterial lengths utilized in this study. distance between the phenotypic and objective curves, which Artery Unit pli,inf pli,sup Reference is computed as follows. The measured waveform fkM (t) is not available in analytic form, but rather as a sequence of values Ascending (i = 2) [cm] 4 5 Goldman and Schafer, 2011 M = f M (t ) at discrete time instants t ∈ [0, T ] as provided fk,s k s s ck Aortic arc (i = 3) [cm] 3.85 5.9 Boufi et al., 2017 by the accelerometer. Since, by construction, all EA-generated Thoracic (i = 4) [cm] 12.9 15.7 Redheuil et al., 2011 curves are defined over the same time interval [0, Tck ] as the Abdominal (i = 5) [cm] 13 16 Drake et al., 2009; Goldman experimentally-measured objective function fkM , the values of and Schafer, 2011 the EA curves at the time instants ts can be easily calculated and Iliac (i = 6) [cm] 3.7 7.5 Bergman, 2007 will be denoted by fk,s P . Then, the Euclidean distance between the Carotid (i = 14) [cm] 20 24.4 Choudhry et al., 2016 functions fkM (t) and fkP (t) is computed as the Euclidean distance M ) and (t , f P ) as d between their discrete versions (ts , fk,s s k,s s 176 cm). In so doing, we leverage the available information on a X M − fP 2 specific subject to obtain an a priori estimate for G y and narrow dk = fk,s k,s for k = 1, . . . , Nc . (13) s its range for the EA search. All baseline values are the same as those reported in Guidoboni et al. (2019) except for Ts , which has The best 100 curves according to the fitness ranking are selected been estimated from the ECG recording as the average of the time as parents for offspring generation. interval between each R peak and the end of the following T wave. STEP 4: OFFSPRING GENERATION. The offspring genotypes This procedure for determining the initial population is are produced as follows: slightly modified when it is applied to each of the remaining selected curves fkM , for k = 2, . . . , Nc . Specifically, g lv , g rv , g vas , and g y are selected as in Equations (10–12, 16), whereas g r and r gir = gi,p r + U (gi,inf r − gi,p r ) , (min(gi,inf r , gi−1 r ) − gi,p ) g l are selected via a uniform distribution within a range of ±3% gir ∈ G r i = {3, 4, 5, 6}, (14) of the fittest genotypes obtained upon the convergence of the EA applied to f1M . Ultimately, this procedure gives the initial gil = gi,p l l + U (gi,inf l − gi,p l ), (gi,sup l − gi,p ) Frontiers in Physiology | www.frontiersin.org 7 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment gil ∈ G l i ∈ I , (15) STEP 6: SOLUTION FEATURES. Once the EA has reached ( y y U (ḡn , α y ḡn ) for α y >1 convergence, the representative features of the EA solutions y are computed as the average over the three best-ranked gn = y y U (α y ḡn , ḡn ) for α y < 1 curves, namely y ∈ Gy, n ∈ N gn (16) 3 3 g lv = gplv + N 0, σ lv P 1 X P ,b P 1 X P ,b AJ,k = AJ,k AK,k = AK,k (20) 3 3 g lv ∈ G lv (17) b=1 b=1 3 3 1 1 P P g rv = gprv + N 0, σ rv P ,b P ,b X X T J,k = TJ,k T K.,k = TK,k . (21) 3 3 rv rv b=1 b=1 g ∈G (18) g vas = gpvas + N 0, σ vas This holds for each EA run pertaining to the k = 1, . . . , Nc g vas ∈ G vas (19) selected curves from the measured BCG. STEP 7: A POSTERIORI BLOOD PRESSURE ESTIMATION. where the subscript p indicates the genetic code of the parent. Ultimately, the EA provides a personalized estimate for the y We note that the interval of variation for gn , with n ∈ N , is cardiovascular parameters summarized in Table 2 which, via the assumed to be the same for parents and offsprings. All offsprings solution of the physiology-based model described in section 2.2, undergo the physiological check described in STEP 2. A total yields an estimate of the distributions of blood volumes and of λ = 3 physiologically-acceptable offsprings are produced blood pressures within the cardiovascular system of a specific by each parent. Finally, the physiologically-acceptable offpsrings subject. To evaluate the reliability of these estimates, we compare and their parents are ranked according to their similarity to the the blood pressure values measured with a cuff directly on the objective curve fkM under consideration and the fittest λ = subject, as described in section 2.1, with the blood pressure values 100 genotypes are selected to form the generation advancing predicted by the model with the personalized parameters. The in the evolution. Utilizing this procedure, we ensure that (i) cuff measures the pressure at the level of the brachial artery the population size is kept constant at M = 300 through the in the arm which, however, is not explicitly included among generations, and that (ii) only the fittest individuals advance from the major arteries of our cardiovascular model (as shown in one generation to the next. section 2.2). To address this issue, we leverage the results of the STEP 5: CONVERGENCE CHECK. The J peak and the K valley Anglo-Cardiff Collaborative Trial, which included approximately are among the most important traits characterizing the BCG 12, 000 individuals across East Anglia and Wales in the United waveform; they are detectable as the most prominent maximum Kingdom (McEniery et al., 2008). The study provides specific and minimum following the R peak in the ECG (Starr and relationships that can be used to estimate the brachial pressure Noordergraaf, 1967). Let us denote by AM M M J,k (resp. AK,k ) and TJ,k from the central aortic pressure. Interestingly, our cardiovascular M (resp. TK,k ) the magnitude of the J peak (resp. K valley) and models provide the central aortic pressure directly as the blood its timing with respect to the preceding R peak calculated for pressure in the ascending aorta (i = 2). In McEniery et al. each of the k selected objective curves, with k = 1, . . . , Nc . (2008), differences in the diastolic values of the central aortic and The location of the J peak and the K valley, along with their brachial pressures were found to be negligible, thereby suggesting magnitudes and timings with respect to the R peaks in the ECG, to assume the two values to be the same. Conversely, the systolic are illustrated in Figure 3. The EA convergence is assessed by brachial values were found to be higher than the corresponding evaluating whether there is an offspring satisfying the following central aortic values, with specific increments and intervals of two criteria: variability provided as a function of age and gender, as shown in Figure 1 of McEniery et al. (2008). Since in this study, we are 1. The predicted J-K magnitude and timings must be within a considering three male subjects in the range of 20–29 years of age, 5% range when compared to those computed for the objective we adopt an increment of 20 mmHg with an interval of variability curve under consideration; of the ± 10 mmHg. 2. The offspring must be within the top Nfit = 20 curves in the fitness ranking (as shown in Step 3). We note that, by requiring the fitness ranking of the offspring 3. RESULTS to be high enough (second criterium), we aim at finding a solution that is close in an average sense to the measured We begin by comparing the BCG curves measured curve while, simultaneously, maximizing the similarity with experimentally with those predicted by the EA algorithm the J-K features (first criterium). If such an offspring exists, (section 3.1). Next, we examine the EA performance in terms convergence is reached, otherwise, the EA continues to the of estimating various parameters in the cardiovascular model next generation. (section 3.2). Finally, the values of PPs estimated by the algorithm A maximum limit of 30 generations has been set as a stopping are compared with those measured with the cuff placed on the criterium in case convergence is not achieved. arm of each subject (section 3.3). Frontiers in Physiology | www.frontiersin.org 8 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment FIGURE 3 | The electrocardiogam (ECG; top) and f M waveform (bottom) acquired synchronously are reported. R peaks in the ECG are marked with red circles and their time location is indicated by red dashed vertical lines. Amplitude and timing for the J peak are reported for the k−th curve (refer to AM M J,k and TJ,k ), whereas the amplitude and timing for the K valley are reported for the following (k + 1)−th curve (refer to AM K,k+1 and T M K,k+1 ). P 3.1. Comparison Between Experimental P 5 M 1 X |AK,k − AK,k | 1(AK , AM K )= × 100 and Predicted BCG Curves 5 |AMK,k | k=1 Figure 4 reports the BCG curves fkM , with k = 1, . . . , 5 measured experimentally (in black) and the corresponding three best- and in the timings of the J and K peaks defined as ranked f P curves predicted via the EA (in colors) obtained for Subject 1. Analogous figures for Subjects 2 and 3 can be 5 P M P 1 X |T J,k − TJ,k | found in the Supplementary Material. Notably, the agreement 1(T J , TJM ) = M| × 100 , between the measured and predicted curves in the systolic part 5 |TJ,k k=1 of the cardiac cycle is quite satisfactory, with a clearly detectable 5 P M similarity in terms of J-K features. During diastole, though, the P 1 X |T K,k − TK,k | 1(T K , TKM ) = M| × 100 . predicted curves are much flatter than the measured curves, 5 |TK,k k=1 capturing only loosely the peaks and valleys that are exhibited experimentally. This result is not unexpected, since the diastolic Remarkably, the mean percent errors in the timings are features of BCG are known to be more challenging to capture approximately one order of magnitude lower than the percent both experimentally and theoretically (Starr and Noordergraaf, errors in the amplitudes. 1967; Guidoboni et al., 2019). A quantitative comparison between the J-K features in the 3.2. Estimation of Cardiovascular experimental and predicted curves for each of the three subjects Parameters via the EA included in the study is summarized in Table 6 by means of We recall that the main output of the EA algorithm is average percent errors in the amplitudes of the J and K peaks the personalized estimate of the physiological and anatomical defined as parameters in G for a given subject. Detailed results are reported in Figure 5 in the case of Subject 1. Analogous figures for Subjects 2 and 3 can be found in the Supplementary Material. For each parameter, we consider the three best-ranked curves 5 P M and we report the mean as the bar height, along with the P 1 X |AJ,k − AJ,k | maximum and the minimum values as black brackets. The 1(AJ , AM J )= × 100 , 5 k=1 |AM J,k | results indicate that the proposed EA algorithm is capable of Frontiers in Physiology | www.frontiersin.org 9 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment FIGURE 4 | Comparison among the BCG curves fkM , with k = 1, . . . , 5 measured experimentally (in black) and the corresponding three best-ranked curves computed via the EA (in colors) for Subject 1. estimating in a consistent manner, over the five selected objective TABLE 6 | Quantitative comparison between the J-K features in the curves, all the parameters characterizing the left ventricle (ELS, experimentally-measured (M) and EA-predicted (P) curves. ELD, ULO, qL , and Ts ) and the arterial Young modulus (E). Features Subject 1 (%) Subject 2 (%) Subject 3 (%) Some of the parameters characterizing the right ventricle can P also be estimated quite consistently (ERD, URO), while others 1(AJ , AM J ) 12.9 15.2 21.1 P show marked differences among the results obtained for the five 1(AK , AM K ) 17.3 11.4 29.6 P curves (qR , ERS). Similar marked differences are displayed by 1(T J , TJM ) 1.7 1.2 1.5 the estimates for peripheral vascular resistance (R7 ). The values P 1(T K , TKM ) 1.5 2.5 3.3 of some estimated parameters (Ts , ULO, E, ERD, URO) turn Percent errors (1) in the amplitudes (A) and the timings (T) are computed for each subject. out to be close to the baseline values reported in Table 3, while others deviate markedly. Interestingly, though, the EA estimates preserve relationships between relative parameter values without explicitly enforcing them, such as the facts that (i) ULO is larger than URO, implying that the capacity for pressure build-up in the curves for the same subject. Since a ground truth for the the left ventricle is larger than that in the right ventricle; and (ii) estimated parameters is not available, we utilized the width of the ELS (resp. ERS) is larger than ELD (resp. ERD), implying that the interval of the estimated parameters as an indicator of the EA end-systolic elastance is larger than the end-diastolic elastance potential for parameter estimation. More precisely, we utilized in both the left and right ventricles. Finally, the estimated values bold fonts to indicate in Table 7 those results for which the for radii, lengths, and locations of the arterial segments are also semidistance between the maximum and minimum is less than reported in Figure 3. They exhibit small differences as a result of 1/3 of the estimated mean value. Interestingly, the estimates the constraints imposed on the genotype generation. of ELD, ULO, Ts , and URO satisfy this criterion for all three A quantitative comparison of the cardiovascular parameters subjects, whereas the estimates of ELS, qL , ERS, qR , and E were estimated via the proposed EA method for the three subjects satisfactory only for two out of three subjects. The value of the included in the study is summarized in Table 7. For each peripheral resistance R7 resulted to be poorly estimated in all estimated parameter, we report the mean value calculated over subjects. It is also worth noticing that the results for all three five objective curves, along with the minimum and maximum subjects confirm that the ULO is estimated to be larger than values (annotated in italics in parenthesis) obtained over all URO, and that ELS (resp. ERS) is estimated to be larger than Frontiers in Physiology | www.frontiersin.org 10 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment FIGURE 5 | Summary of the physiological and anatomical parameters estimated by the EA for each of the k = 1, . . . , 5 selected BCG objective curves for Subject 1. ELD (resp. ERD), thereby supporting the physiological relevance 3 can be found in the Supplementary Material. Figure 6 (Left) of the findings. shows the pressure waveforms in the ascending aorta (i = 2) predicted by the cardiovascular model for the three best-ranked 3.3. Central Aortic Pressure and Brachial curves obtained by the EA performed on f1M . Following McEniery Pressure et al. (2008), we apply a 20 mmHg increment (vertical yellow To further verify the capability of the EA algorithm to yield segments) to the predicted systolic value of the central aortic physiologically-meaningful solutions, we compare the brachial pressure to estimate the systolic value of the brachial pressure. pressure measured experimentally with the pressure predicted Differences in the diastolic values of the central aortic and by the cardiovascular model equipped with the personalized brachial pressure are neglected. Figure 6 (Right) reports the parameters provided by the EA search. The results obtained results for the brachial pressure predicted by the cardiovascular for Subject 1 are reported in Figure 6, where the horizontal model for Subject 1 with the personalized model parameters lines indicate the mean (solid line), maximum, and minimum yielded by the EA search performed on each of the k = 1, . . . , 5 values (dashed lines) of the six systolic and diastolic pressure objective curves fkM . The height of the colored bars represents measurements acquired with a cuff placed on the arm of the the mean value over the three best-ranked curves obtained for a subject (refer to section 2.1). Analogous figures for Subjects 2 and given k, whereas the black brackets indicate the maximum and Frontiers in Physiology | www.frontiersin.org 11 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment TABLE 7 | Summary of EA estimated parameters for the three subjects involved in the study. Parameter Unit Subject 1 Subject 2 Subject 3 ELS mmHg cm−3 0.37 (0.32,0.43) 0.37 (0.28,0.48) 0.78 (0.50,1.50) ELD mmHg cm−3 0.08 (0.07,0.10) 0.09 (0.08,0.10) 0.06 (0.04,0.08) ULO mmHg 66.92 (62.22,71.84) 80.22 (74.21,87.31) 74.23 (58.11,90.89) qL s−1 13.54 (12.25,16.10) 12.85 (10.48,16.02) 15.76 (3.51,26.08) Ts s 0.37 (0.35,0.38) 0.37 (0.36,0.38) 0.37 (0.36,0.39) ERS mmHg cm−3 0.42 (0.32,0.60) 0.37 (0.30,0.52) 0.37 (0.26,0.42) ERD mmHg cm−3 0.019 (0.011,0.024) 0.018 (0.012,0.023) 0.020 (0.010,0.037) URO mmHg 20.15 (11.76,24.90) 21.73 (16.58,26.16) 23.95 (18.28,30.67) qR s−1 9.98 (6.62,15.66) 8.36 (6.48,10.42) 9.71 (7.84,11.34) 3 E 10 mmHg 2.60 (2.30,2.97) 4.04 (3.45,4.88) 4.60 (3.33,8.22) R7 mmHg cm−3 s 0.14 (0.09,0.19) 0.19 (0.14,0.27) 0.11 (0.06,0.16) Mean values are reported along with the minimum and maximum values (in italics, in parenthesis). Bold fonts are used to indicate the instances for which the semidistance between maximum and minimum is less than 1/3 than the estimated mean value. minimum values. The pulse pressure, defined as the difference (refer to Figure 4 and Table 6). By doing so, we are able to between systolic and diastolic values, is highlighted with solid obtain personalized estimates of cardiovascular parameters (refer colors in the bars. to Figure 5 and Table 7) and of variables of physiological interest, A comparison between the values of PP measured such as the central aortic and brachial pressures (see Figure 6 and experimentally and predicted by the proposed EA method Table 8). for each subject of three subjects is summarized in Table 8. In the current implementation of the algorithm, we opted for The measured values correspond to the average PP values selecting Nc = 5 consecutive experimental BCG curves rather over a total of six measurements obtained with the cuff, as than a single template representing the whole data acquisition. described in section 2.1. The values reported in italics in This choice is motivated by the fact that amplitude, timing, and parenthesis indicate the minimum and maximum values length of each BCG curve embody the action and function of in the single measurements. The EA-predicted values are the ventricles and the vasculature during a single heart-beat and obtained through the following steps: (i) a central aortic may as well vary in the next. Thus, we do not expect the EA pressure waveform (refer to Figure 6) is obtained via the estimates for the personalized cardiovascular parameters to be cardiovascular model (refer to section 2.2) with the set of model the same from beat-to-beat, even for healthy subjects such as parameters corresponding to the three best-ranked curves; those considered in this study, but rather to be in the same (ii) a factor of 20 ± 10 mmHg is added to the systolic value ballpark, as shown in Figure 5. A limitation of the current study of the simulated central aortic pressure to obtain the systolic is that the values of the estimated parameters are not directly brachial pressure, as suggested by the population-based study comparable to independent measurements. Radii, length, and of McEniery et al. (2008) (refer to section 2.1); (iii) the PP locations of the main arteries, G r , G l , and G y could be acquired, is calculated as the difference between the estimated systolic for example, using Doppler imaging. Such information could brachial pressure and the diastolic blood pressure; (iv) the either be used as a posterior verification of the EA predictions or PP values are averaged over the three best-ranked curves for as a priori knowledge that would narrow the EA search range for each of the five objective curves for each subject, with the that specific subject. The measurement of other parameters, such overall minimum and maximum values reported in italics in as the ventricular elastances, requires invasive techniques based parenthesis. Since a ground truth for the central aortic pressure on catheterization. The capability of the proposed approach is not available for this study, the fact that the predicted PP values to yield physiologically-meaningful solutions is confirmed by are within the measured intervals for all three subjects is very the good agreement between the PP values predicted by the promising and provides supportive indirect evidence that the cardiovascular model personalized via the EA method and the parameters estimated via the proposed EA method actually bear experimental measurements (refer to Figure 6 and Table 8). This physiological relevance. result is particularly encouraging, considering that the blood pressure values were not included in any of the feature sets 4. DISCUSSIONS utilized in the EA for the physiological check (refer to Step 2, section 2.3), fitness evaluation (refer to Step 3, section 2.3), and The novelty of the approach proposed in the present study convergence (refer to Step 5, section 2.3). This finding suggests consists in leveraging a physiology-based mathematical model that our approach could be used to obtain BCG-based cuffless to incorporate substantial domain knowledge in an EA whose blood pressure measurements, along with noninvasive estimates objective is to attain optimal fitness between model-predicted of central aortic pressure and valuable parameters describing and experimentally-measured BCG curves on a given subject cardiovascular function. Frontiers in Physiology | www.frontiersin.org 12 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment FIGURE 6 | Left: Comparison between the central aortic pressure corresponding to the three best-ranked curves selected by the EA search performed on f1M for Subject 1 (yellow curves) and the blood pressure measured at the arm with a cuff of the subject (horizontal black lines). The 20 mmHg increment applied to the systolic value of predicted central aortic pressure is also indicated (vertical yellow segments). The mean (solid black lines) maximum and minimum value (dashed black lines) of the repeated blood pressure measurements are reported. Right: Comparison between the experimentally measured brachial pulse pressure (PP) (horizontal blue lines) and the brachial pressure predicted by the EA for each of the k = 1, . . . , 5 objective curves for Subject 1. The PP is indicated with solid colors. Maximum and minimum values obtained for the three best-ranked curves for each fkM , with k = 1, . . . , Nc , are reported in black brackets. In the long term, this study aims at contributing to the quest TABLE 8 | Comparison between the pulse pressure (PP) values measured via a for noninvasive techniques capable of providing meaningful cuff placed in the arm and the PP values predicted as a result of the cardiovascular parameters estimated via the proposed EA method. insights into the thermodynamic efficiency of cardiac function. Without direct left ventricular inductance catheters, clinicians Pulse Pressure [mmHg] Subject 1 Subject 2 Subject 3 must rely on indirect estimations from right heart catheters or algorithms from echocardiograms, each affected by risks and Measured 57.5 40.5 57.3 limitations (Ikonomidis et al., 2019). By providing noninvasive (53,62) (35,47) (52,60) estimates of left-ventricular end-systolic elastance and central Predicted 61.1 41.6 44.7 aortic pressure based on a mechanistic interpretation of the (51.1,71.1) (31.6,51.6) (34.7,54.7) BCG signal, this study could provide clinicians with a rapid Variability intervals are reported in italics in parenthesis. and insightful assessment of cardiac function that could be used at the bedside of the critically ill patient and offer practical solutions for outpatient monitoring. To get a sense of how these results could be used in practice, let us look at the parameters currently conducting studies on human subjects and on swine to estimated for Subjects 1 and 2 in Table 7. Based on the EA- provide further data supporting the theoretical findings reported guided interpretation of the BCG signal, the pressure build-up in the article, hopefully, bring us closer to making this vision capacity in the left ventricle (ULO) for Subject 2 is approximately a reality. 20% higher than in Subject 1, while being very close for the Limitations from both the experimental and modeling right ventricle (URO). This difference does not constitute a viewpoints should be considered when evaluating the findings problem per se; rather, it shows how the BCG could be used of our study. The BCG sensing modality utilized in this study to establish a cardiovascular baseline for each individual. In the is a suspended bed equipped with an accelerometer. While case of outpatient monitoring, longitudinal measurements over providing a signal that is very close to the true acceleration of the course of months could help detect a deterioration in left- the center of mass of the human body (Starr and Noordergraaf, ventricular function by, for example, providing a quantitative 1967), this sensing modality is primarily used in research trajectory of decreasing ULO values. In the case of critically laboratories and is not amenable to clinical or in-home use. ill patients, frequent monitoring (possibly continuous) may be Despite differing in shape among sensing modalities, all BCG advisable in order to enable early detection of cardiogenic shock. waveforms exhibit a major peak (i.e., J peak) and a major Similarly, this method could be used to track changes in left- valley (i.e., K valley) (Giovangrandi et al., 2011). Thus, by using ventricular end-systolic and end-diastolic elastances (ELS, ELD), only the J-K amplitudes and timings in the convergence of the whose changes are indicative of heart failure with reduced and EA algorithm, the approach described in this study could be preserved EF, respectively (Guidoboni et al., 2019). Our group is extended to other BCG technologies. An ongoing study in the Frontiers in Physiology | www.frontiersin.org 13 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment Surgical Intensive Care Unit (MU Health Care, Columbia, MO) a cuff shows that the proposed approach is physiologically has recently shown that measures of timing between ECG and meaningful and may provide theoretical support to the BCG signals acquired on critically ill patients by means of a further development of cuffless methods for blood pressure three-axis accelerometer positioned under the head pillow are measurements (Solà and Delgado-Gonzalo, 2019; Le et al., feasible and reproducible (Zaid et al., 2021), thereby showing 2020; Pandit et al., 2020). Investigations evaluating the good potential for applications of the proposed methodology applicability of the proposed approach to situations where beyond a laboratory setting. the data acquisition is not as controlled as in a laboratory Additionally, our results show that the agreement between setting are currently under-way. Preliminary results obtained model-predicted and experimentally measured BCG curves is when monitoring critically ill patients hospitalized in the better in the systolic part than in the diastolic part of the Surgical Intensive Care Unit (University Hospital, MU cardiac cycle (refer to Figure 4). It is known that the BCG signal Health Care System) suggest that measurements of BCG is stronger during systole when the ventricular contractions amplitudes and timings are feasible and reproducible (Zaid occur and the blood from the left ventricle is ejected and et al., 2021), thereby yielding promise for future extension of channeled through the aorta (Starr and Noordergraaf, 1967; this study. Kim et al., 2016). Thus, the experimental measurements are much more reliable during systole than diastole. Furthermore, DATA AVAILABILITY STATEMENT the physiological-based cardiovascular model for BCG prediction used in this study is capable of simulating the systolic peak The original contributions presented in the study are included and valleys of the BCG waveform with much greater accuracy in the article/Supplementary Material, further inquiries can be than those in the diastole (Guidoboni et al., 2019). Due to directed to the corresponding author. these experimental and theoretical limitations, we based our convergence criteria on systolic features of the BCG waveform. ETHICS STATEMENT In future studies, with the advances of BCG technologies and physiological understanding of the BCG waveform, these features The studies involving human participants were reviewed and could be extended to include also the diastolic part of the cardiac approved by Institutional Review Board, Office of Research cycle. An aspect that could be considered in evaluating the and Economic Development, University of Missouri. The performance of the proposed EA method is a different choice patients/participants provided their written informed consent to for Nc representing the number of consecutive BCG curves to participate in this study. be selected as objective curves. The choice of Nc = 5 adopted in this study is motivated by the need of considering multiple curves AUTHOR CONTRIBUTIONS while maintaining the overall computational load affordable. Optimal choices for Nc may be explored in conjunction with NM, LS, and MZ contributed to implementing the algorithm the effect of breathing, which may affect the BCG curves over and running the numerical simulations. GG, NM, and longer intervals. MZ wrote the first draft of the manuscript. All authors contributed to the conception and design of the study, 5. CONCLUSIONS contributed to manuscript revision, read, and approved the submitted version. This study presented a novel combination of a physiology- based mathematical model and an evolutionary algorithm to ACKNOWLEDGMENTS obtain personalized estimates of cardiovascular parameters and variables of physiological interest, such as blood pressure, The authors acknowledge support from the Center of Eldercare with the goal of developing quantitative tools for noninvasive and Rehabilitation Technology and the University of Missouri. cardiovascular evaluations based on BCG sensing. The approach proved capable of estimating many ventricular and arterial SUPPLEMENTARY MATERIAL parameters with consistency when five consecutive BCG curves were selected for the subjects considered in this study. The Supplementary Material for this article can be found Furthermore, the good agreement between the blood pressure online at: https://www.frontiersin.org/articles/10.3389/fphys. estimated with the model and measured experimentally with 2021.739035/full#supplementary-material REFERENCES Bergman, R. A. (2007). Illustrated Encyclopedia of Human Anatomic Variation. Available online at: http://www. anatomyatlases. ADInstruments (2014). Powerlab/35 series. Available online at: https://www. org/AnatomicVariants/MuscularSystem/Text/P/08Pectoralis. shtml. adinstruments.com/products/powerlab/35-and-26-series Boufi, M., Guivier-Curien, C., Loundou, A., Deplano, V., Boiron, O., Alametsä, J., Viik, J., Alakare, J., Vrri, A., and Palomki, A. (2008). Chaumoitre, K., et al. (2017). Morphological analysis of healthy aortic Ballistocardiography in sitting and horizontal positions. Physiol. Meas. arch. Eur. J. Vasc. Endovasc. Surg. 53, 663–670. doi: 10.1016/j.ejvs.2017. 29, 1071. doi: 10.1088/0967-3334/29/9/006 02.023 Frontiers in Physiology | www.frontiersin.org 14 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment Chang, M. C., Martin, R. S., Scherer, L. A., and Meredith, J. W. (2002). Improving Kim, C.-S., Ober, S. L., McMurtry, M. S., Finegan, B. A., Inan, O. T., Mukkamala, ventricular-arterial coupling during resuscitation from shock: effects on R., et al. (2016). Ballistocardiogram: Mechanism and potential for unobtrusive cardiovascular function and systemic perfusion. J. Trauma Acute Care Surg. cardiovascular health monitoring. Sci. Rep. 6:31297. doi: 10.1038/srep31297 53, 679–685. doi: 10.1097/00005373-200210000-00010 Kionix, I. (2014). Kxr94-2283. Available online at: https://www.kionix.com/ Chen, W., Zhu, X., Nemoto, T., Kitamura, K., Sugitani, K., and Wei, D. (2008). Krejza, J., Arkuszewski, M., Kasner, S. E., Weigele, J., Ustymowicz, A., Unconstrained monitoring of long-term heart and breath rates during sleep. Hurst, R. W., et al. (2006). Carotid artery diameter in men and Physiol. Meas. 29, N1. doi: 10.1088/0967-3334/29/2/N01 women and the relation to body and neck size. Stroke 37, 1103–1105. Choudhry, F. A., Grantham, J. T., Rai, A. T., and Hogg, J. P. (2016). doi: 10.1161/01.STR.0000206440.48756.f7 Vascular geometry of the extracranial carotid arteries: an analysis Le, T., Ellington, F., Lee, T.-Y., Vo, K., Khine, M., Krishnan, S., et al. (2020). of length, diameter, and tortuosity. J. Neurointerv. Surg. 8, 536–540. Continuous non-invasive blood pressure monitoring: a methodological doi: 10.1136/neurintsurg-2015-011671 review on measurement techniques. IEEE Access. 8, 212478–212498. Drake, R., Vogl, A. W., and Mitchell, A. W. (2009). Gray’s Anatomy for Students E- doi: 10.1109/ACCESS.2020.3040257 Book. Elsevier Health Sciences. Available online at: https://www.elsevier.com/ McEniery, C. M., Yasmin, n., McDonnell, B., Munnery, M., Wallace, S. M., Rowe, books/grays-anatomy-for-students/drake/978-0-323-39304-1 C. V., et al. (2008). Central pressure: variability and impact of cardiovascular Enayati, M. (2019). Data-Driven Methods for Analyzing Ballistocardiograms in risk factors: the anglo-cardiff collaborative trial ii. Hypertension 51, 1476–1482. Longitudinal Cardiovascular Monitoring (Dissertation). University of Missouri, doi: 10.1161/HYPERTENSIONAHA.107.105445 USA. Paalasmaa, J., Toivonen, H., and Partinen, M. (2014). Adaptive heartbeat modeling Etemadi, M., Inan, O. T., Giovangrandi, L., and Kovacs, G. T. (2011). Rapid for beat-to-beat heart rate measurement in ballistocardiograms. IEEE J. Biomed. assessment of cardiac contractility on a home bathroom scale. IEEE Trans. Inf. Health Inform. 19, 1945–1952. doi: 10.1109/JBHI.2014.2314144 Technol. Biomed. 15, 864–869. doi: 10.1109/TITB.2011.2161998 Pan, J., and Tompkins, W. J. (1985). A real-time qrs detection algorithm. IEEE Giovangrandi, L., Inan, O. T., Wiard, R. M., Etemadi, M., and Kovacs, G. T. (2011). Trans. Biomed. Eng. 32, 230–236. doi: 10.1109/TBME.1985.325532 Ballistocardiography a method worth revisiting. Annu. Int. Conf. IEEE Eng. Pandit, J. A., Lores, E., and Batlle, D. (2020). Cuffless blood pressure Med. Biol. Soc. 2011, 4279–4282. doi: 10.1109/IEMBS.2011.6091062 monitoring: promises and challenges. Clin. J. Am. Soc. Nephrol. 15, 1531–1538. Goldman, L., and Schafer, A. I. (2011). Goldman’s cecil medicine E-book. doi: 10.2215/CJN.03680320 Elsevier Health Sciences. Available online at: https://www.elsevier.com/books/ Ranganathan, S., Nakai, K., and Schonbach, C. (2018). Encyclopedia of goldman-cecil-medicine-2-volume-set/goldman/978-0-323-53266-2 Bioinformatics and Computational Biology: ABC of Bioinformatics. Elsevier. Guidoboni, G. (2020). Model-Based Sensor Technology for Detection of Available online at: https://www.elsevier.com/books/encyclopedia-of- Cardiovascular Status. Assignee: The Curators of the University of bioinformatics-and-computationalbiology/ranganathan/978-0-12-811414-8 Missouri. PCT Patent Application No.: PCT/US2019/052738. Published Redheuil, A., Yu, W.-C., Mousseaux, E., Harouni, A. A., Kachenoura, N., Wu, C. as WO2020068859A1. O., et al. (2011). Age-related changes in aortic arch geometry: relationship with Guidoboni, G., Sala, L., Enayati, M., Sacco, R., Szopos, M., Keller, J. M., et al. (2019). proximal aortic function and left ventricular mass and remodeling. J. Am. Coll. Cardiovascular function and ballistocardiogram: a relationship interpreted Cardiol. 58, 1262–1270. doi: 10.1016/j.jacc.2011.06.012 via mathematical modeling. IEEE Trans. Biomed. Eng. 66, 2906–2917. Sacco, R., Guidoboni, G., and Mauri, A. G. (2019). A Comprehensive Physically doi: 10.1109/TBME.2019.2897952 Based Approach to Modeling in Bioengineering and Life Sciences. Academic Heise, D., Rosales, L., Skubic, M., and Devaney, M. (2011). Refinement and Press. Available online at: https://www.elsevier.com/books/a-comprehensive- evaluation of a hydraulic bed sensor. Annu. Int. Conf. IEEE Eng. Med. Biol. Soc. physically-based-approach-to-modeling-inbioengineering-and-life-sciences/ 2011, 4356–4360. doi: 10.1109/IEMBS.2011.6091081 sacco/978-0-12-812518-2 Helfand, M., Christensen, V., and Anderson, J. (2016). “Technology assessment: Satu, R., and Jukka, L. (2012). Film-type sensor materials PVDF and EMFi in early sense for monitoring vital signs in hospitalized patients,” in VA Evidence measurement of cardiorespiratory signals. a review. IEEE Sens. J. 12, 439–446. Synthesis Program Evidence Briefs [Internet]. Washington, DC: Department of doi: 10.1109/JSEN.2010.2089510 Veterans Affairs (US). Shin, J., Choi, B., Lim, Y., Jeong, D., and Park, K. (2008). Automatic Holcroft, J. W., John, A. T., and Sena, M. J. (2006). Management of shock. ACS ballistocardiogram (BCG) beat detection using a template matching approach. Surg. Princ. Pract. 3, 1–39. Annu. Int. Conf. IEEE Eng. Med. Biol. Soc. 2008, 1144–1146. Huffaker, M., Carchia, M., Harris, B., Kethman, W., Murphy, T., Sakarovitch, C., Solà, J., and Delgado-Gonzalo, R. (2019). The Handbook of Cuffless Blood Pressure et al. (2018). Passive nocturnal physiologic monitoring enables early detection Monitoring. Springer. Available online at: https://link.springer.com/book/10. of exacerbations in children with asthma. A proof-of-concept study. Am. J. 1007/978-3-030-24701-0 Respir. Crit. Care Med. 198, 320–328. doi: 10.1164/rccm.201712-2606OC Soucier, R., Miller, P., Ingrassia, J., Riello, R., Desai, N., and Ahmad, T. (2018). Ikonomidis, I., Aboyans, V., Blacher, J., Brodmann, M., Brutsaert, D. L., Chirinos, Essential elements of early post discharge care of patients with heart failure. J. A., et al. (2019). The role of ventricular-arterial coupling in cardiac Curr. Heart Fail. Rep. 15, 181–190. doi: 10.1007/s11897-018-0393-9 disease and heart failure: assessment, clinical implications and therapeutic Starr, I., and Noordergraaf, A. (1967). Ballistocardiography in Cardiovascular interventions. a consensus document of the european society of cardiology Research: Physical aspects of the Circulation in Health and Disease. working group on aorta peripheral vascular diseases, european association of Lippincott. Available online at: https://books.google.com/books/about/ cardiovascular imaging, and heart failure association. Eur. J. Heart Fail. 21, Ballistocardiography_in_Cardiovascular_R.html?id=_iBsAAAAMAAJ 402–424. doi: 10.1002/ejhf.1436 Su, B. Y., Enayati, M., Ho, K., Skubic, M., Despins, L., Keller, J., et al. (2018). Inan et al., O. (2009). Robust ballistocardiogram acquisition for home monitoring. Monitoring the relative blood pressure using a hydraulic bed sensor system. Physiol. Meas. 30, 169. doi: 10.1088/0967-3334/30/2/005 IEEE Trans. Biomed. Eng. 66, 740–748. doi: 10.1109/TBME.2018.2855639 Inan, O. T., Migeotte, P.-F., Park, K.-S., Etemadi, M., Tavakolian, K., Vincent, J.-L., and De Backer, D. (2013). Circulatory shock. N. Engl. J. Med. 369, Casanella, R., et al. (2014). Ballistocardiography and seismocardiography: a 1726–1734. doi: 10.1056/NEJMra1208943 review of recent advances. IEEE J. Biomed. Health Inform. 19, 1414–1427. Wasfy, J., Strom, J., Waldo, S., O’Brien, C., Wimmer, N., Zai, A., et al. (2014). doi: 10.1109/JBHI.2014.2361732 Clinical preventability of 30-day readmission after percutaneous coronary Joh, J. H., Ahn, H.-J., and Park, H.-C. (2013). Reference diameters of the abdominal intervention. J. Am. Heart Assoc. 3:e001290. doi: 10.1161/JAHA.114.001290 aorta and iliac arteries in the korean population. Yonsei Med. J. 54, 48–54. Wolak, A., Gransar, H., Thomson, L. E., Friedman, J. D., Hachamovitch, doi: 10.3349/ymj.2013.54.1.48 R., Gutstein, A., et al. (2008). Aortic size assessment by noncontrast Katz, Y., Karasik, R., and Shinar, Z. (2016). “Contact-free piezo electric sensor cardiac computed tomography: normal limits by age, gender, and body used for real-time analysis of inter beat interval series,” in 2016 Computing in surface area. JACC Cardiovas. Imaging 1, 200–209. doi: 10.1016/j.jcmg.2007. Cardiology Conference (CinC) (Vancouver, BC: IEEE), 769–772. 11.005 Frontiers in Physiology | www.frontiersin.org 15 January 2022 | Volume 12 | Article 739035 Marazzi et al. Modeling for Personalized Cardiovascular Assessment World Health O.rganization (2021). Cardiovascular Diseases (CVDs). Available The remaining authors declare that the research was conducted in the absence of online at: https://www.who.int/news-room/fact-sheets/detail/cardiovascular- any commercial or financial relationships that could be construed as a potential diseases-(cvds). conflict of interest. Young, S., Gillon, W., Rifredi, R., and Krein, W. (2008). Method and apparatus for monitoring vital signs remotely. US Patent App. 11/849,051. Zaid, M., Ahmad, S., Suliman, A., Camazine, M., Weber, I., Sheppard, J., et al. Publisher’s Note: All claims expressed in this article are solely those of the authors (2021). Noninvasive cardiovascular monitoring based on electrocardiography and do not necessarily represent those of their affiliated organizations, or those of and ballistocardiography: a feasibility study on patients in the surgical intensive the publisher, the editors and the reviewers. Any product that may be evaluated in care unit. Annu. Int. Conf. IEEE Eng. Medicine Biol. Soc. 2021, 951–954. this article, or claim that may be made by its manufacturer, is not guaranteed or doi: 10.1109/EMBC46164.2021.9629531 endorsed by the publisher. Zimlichman, E., Szyper-Kravitz, M., Shinar, Z., Klap, T., Levkovich, S., Unterman, A., et al. (2012). Early recognition of acutely deteriorating patients in non- intensive care units: Assessment of an innovative monitoring technology. J. Copyright © 2022 Marazzi, Guidoboni, Zaid, Sala, Ahmad, Despins, Popescu, Hosp. Med. 7, 628–633. doi: 10.1002/jhm.1963 Skubic and Keller. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution Conflict of Interest: GG would like to disclose that she received remuneration or reproduction in other forums is permitted, provided the original author(s) from Foresite Healthcare LLC for serving as a consultant. MS also discloses a and the copyright owner(s) are credited and that the original publication in conflict with Foresite Healthcare LLC outside the submitted work and patents this journal is cited, in accordance with accepted academic practice. No use, licensed to Foresite Healthcare LLC. These relationships are pursuant to the distribution or reproduction is permitted which does not comply with these University of Missouri’s policy on outside activities. terms. Frontiers in Physiology | www.frontiersin.org 16 January 2022 | Volume 12 | Article 739035