@prefix vivo: . @prefix edm: . @prefix ns0: . @prefix dcterms: . @prefix skos: . vivo:departmentOrSchool "Applied Science, Faculty of"@en, "Electrical and Computer Engineering, Department of"@en ; edm:dataProvider "DSpace"@en ; ns0:degreeCampus "UBCV"@en ; dcterms:creator "Leung, Mande Tak Man"@en ; dcterms:issued "2011-06-30T20:06:41Z"@en, "2011"@en ; vivo:relatedDegree "Doctor of Philosophy - PhD"@en ; ns0:degreeGrantor "University of British Columbia"@en ; dcterms:description """According to the World Health Organization (WHO), cardiovascular (CV) diseases are the number one cause of death globally, and are projected to remain the single leading cause of death by 2030. Central pulse pressure and several arterial stiffness indices, especially aortic pulse wave velocity (PWV), are among the strongest predictors of CV events, including stroke, myocardial infarction (heart attack) and angina (chest pain). Therefore, non-invasive methods for the assessment of the central pressure and central arterial stiffness are very important for the diagnosis of CV diseases at their early stage of development. In this thesis, non-invasive methods are developed to (i) measure the aortic PWV, and (ii) measure the aortic pressure waveform (APW). In the intensive care unit, cardiac output (CO) is an important measure of the adequacy of circulation among post-operative patients of cardiovascular diseases. Current methods of monitoring use the heart rate, peripheral blood pressure, and urine output as surrogates of CO, but these indices are inadequate. CO can be monitored directly using thermodilution but the procedure is highly invasive. In this thesis, a minimally-invasive method is developed to monitor the CO using the radial artery pressure waveform. The underlying algorithms of the methods developed in this thesis are interrelated. In this thesis, the arterial system is modeled as a single uniform lossless transmission line terminated by a complex load. Each measuring and monitoring task is a problem of identification, simulation or data acquisition. The measurement of the aortic PWV is formulated as an identification problem. The method developed in this thesis was shown to be able to differentiate a group of children with Marfan syndrome from the healthy children. The monitoring of CO is formulated as a problem of identification and simulation. The method developed in this thesis was applied to five post-surgical infants. It showed clinically acceptable agreement with the more established echocardiographic technique. The measurement of APW is a problem of data acquisition. A method was developed to estimate APW from the aortic distension waveform obtained using B-mode ultrasound. This method showed good agreement with carotid artery applanation tonometry when applied to nine healthy children."""@en ; edm:aggregatedCHO "https://circle.library.ubc.ca/rest/handle/2429/35847?expand=metadata"@en ; skos:note "The Measurement and the Monitoring of the Central Hemodynamics using Non-invasive and Minimally-Invasive Data by Mande Tak Man Leung B.Sc., The University of Alberta, 2004 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in The Faculty of Graduate Studies (Electrical and Computer Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) June 2011 c© Mande Tak Man Leung 2011 Abstract According to the World Health Organization (WHO), cardiovascular (CV) diseases are the number one cause of death globally, and are projected to remain the single leading cause of death by 2030. Central pulse pressure and several arterial stiffness indices, especially aortic pulse wave velocity (PWV), are among the strongest predictors of CV events, including stroke, myocardial infarction (heart attack) and angina (chest pain). Therefore, non-invasive methods for the assessment of the central pressure and central arterial stiffness are very important for the diagnosis of CV diseases at their early stage of development. In this thesis, non-invasive methods are developed to (i) measure the aortic PWV, and (ii) measure the aortic pressure waveform (APW). In the intensive care unit, cardiac output (CO) is an important measure of the adequacy of circulation among post-operative patients of cardiovascu- lar diseases. Current methods of monitoring use the heart rate, peripheral blood pressure, and urine output as surrogates of CO, but these indices are inadequate. CO can be monitored directly using thermodilution but the procedure is highly invasive. In this thesis, a minimally-invasive method is developed to monitor the CO using the radial artery pressure waveform. The underlying algorithms of the methods developed in this thesis are interrelated. In this thesis, the arterial system is modeled as a single uniform lossless transmission line terminated by a complex load. Each measuring and monitoring task is a problem of identification, simulation or data acquisition. The measurement of the aortic PWV is formulated as an identification problem. The method developed in this thesis was shown to be able to differentiate a group of children with Marfan syndrome from the healthy children. The monitoring of CO is formulated as a problem of identification and simulation. The method developed in this thesis was applied to five post- surgical infants. It showed clinically acceptable agreement with the more established echocardiographic technique. The measurement of APW is a problem of data acquisition. A method ii Abstract was developed to estimate APW from the aortic distension waveform ob- tained using B-mode ultrasound. This method showed good agreement with carotid artery applanation tonometry when applied to nine healthy children. iii Preface This thesis is based on several manuscripts, resulting from collaboration between multiple researchers. A version of Chapter 2 has been submitted. This paper was co-authored with Prof Guy Dumont, Dr George Sandor and Dr James Potts. The author’s contribution was developing the idea, numerical implementation, processing the raw echocardiographic data and tonometric data, algorithm evaluation, comparing results between the normal subjects and diseased sub- jects ,and preparing the manuscript. Prof Dumont assisted with his insight on the identification technique and optimization algorithm, and edited the manuscript. Dr. Sandor assisted with his insight on the potential clinical use of the technique and processed the echocardiographic data. Dr. Potts provided his expertise on statistics for medical research. A version of Chapter 3 has been submitted. This paper was co-authored with Prof Guy Dumont, Dr George Sandor and Dr James Potts. The au- thor’s contribution was developing the idea, numerical implementation, col- lecting radial artery pressure data in the intensive care unit, developing a software for data extraction in the intensive care unit, processing the raw echocardiographic and pressure data, algorithm evaluation using simulated data, technique evaluation using data collected from human subjects in the intensive care unit, and preparing the manuscript. Prof Dumont assisted with his insight on the identification technique, suggested reviews on other recently developed similar techniques and edited the manuscript. Dr. San- dor assisted with his insight on the potential clinical use of the technique, assisted with subject selections, evaluated the conditions of the subjects prior to and during data collection, commented on the performance of the technique in a clinical setting. Dr. Potts provided his expertise on statistics for medical research. A version of Chapter 4 is under preparation for submission. This paper was co-authored with Prof Guy Dumont, Dr George Sandor and Dr James Potts. The author’s contribution was developing the idea, numerical imple- mentation, processing the raw echocardiographic data and tonometric data, comparing the proposed technique with the more established tonometric iv Preface technique, and preparing the manuscript. Prof Dumont assisted with his insight on the mathematical aspect of the technique. Dr. Sandor assisted with his insight on the potential clinical use of the technique. Dr. Potts provided his expertise on statistics for medical research. A version of Chapter 2 was accepted and presented as • M. Leung, G. A. Dumont, J. E. Potts, M. T. Potts and G. G. S. San- dor, “A model-based method to estimate pulse wave velocity,” 11th World Congress of Medical Physics and Biomedical Engineering, Mu- nich, Germany, 2009 Another version of Chapter 2 was accepted and presented as • M. Leung, G. A. Dumont, G. G. S. Sandor and J. E. Potts, “Esti- mating arterial stiffness using transmission line model,” 28th IEEE Engineering in Medicine and Biology Conference, New York, 2006 A version of Chapter 4 was accepted and presented as • M. Leung, G. A. Dumont, G. G. S. Sandor, M. T. Potts and J. E. Potts, “A novel method to estimate the aortic pressure waveform using B-mode ultrasound images acquired from a suprasternal view,” 30th IEEE Engineering in Medicine and Biology Conference, Vancouver, 2008 The data collection for the studies presented in Chapter 2 and Chap- ter 4 was approved by the University of British Columbia Children’s and Women’s Research Ethics Board. The certificate number is H07-01816. The data collection for the study presented in Chapter 3 was approved by the University of British Columbia Children’s and Women’s Research Ethics Board. The certificate number is H06-038876. v Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Preface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iv Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xi Abbreviations and Symbols . . . . . . . . . . . . . . . . . . . . . xiii Dedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xviii 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2.1 The mechanism of the pressure and flow waves in the arterial system . . . . . . . . . . . . . . . . . . . . . . 3 1.3 Literature review . . . . . . . . . . . . . . . . . . . . . . . . 5 1.3.1 Current non-invasive methods to measure the central pressure waveform . . . . . . . . . . . . . . . . . . . . 5 1.3.2 Current non-invasive methods to measure aortic pulse wave velocity . . . . . . . . . . . . . . . . . . . . . . . 8 1.3.3 Current methods to monitor cardiac output . . . . . 11 1.4 Thesis objectives . . . . . . . . . . . . . . . . . . . . . . . . . 14 2 Model-based Aortic Pulse Wave Velocity1,2,3 . . . . . . . . 17 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 2.2.1 Outline of the model parameter identification algo- rithm . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 vi Table of Contents 2.2.2 Estimations of the line inductance and the terminat- ing resistance . . . . . . . . . . . . . . . . . . . . . . 21 2.2.3 Initial values of the line capacitance, the load capaci- tance and the line length . . . . . . . . . . . . . . . . 21 2.2.4 Estimation of the load resistance . . . . . . . . . . . . 22 2.2.5 Performance index . . . . . . . . . . . . . . . . . . . . 22 2.2.6 Calculation of the model-based pulse wave velocity . 24 2.3 Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 2.3.1 Data acquisition and processing . . . . . . . . . . . . 24 2.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 2.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 2.5.1 Ability to differentiate the subjects with Marfan syn- drome from the healthy subjects . . . . . . . . . . . . 35 2.5.2 Correlations between the pulse wave velocities obtained by different methods . . . . . . . . . . . . . . . . . . 35 2.5.3 Interchangeability between the pulse wave velocities obtained by different methods . . . . . . . . . . . . . 36 2.5.4 Performance of the model in predicting ascending aor- tic flow . . . . . . . . . . . . . . . . . . . . . . . . . . 36 2.5.5 Large artery compliance to small artery compliance ratio . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 2.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3 Monitoring Cardiac Output using Grey-Box Identification on a State-Space Representation of the Transmission Line Model1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 3.2.1 Derivation of the state-space representation of the trans- mission line model . . . . . . . . . . . . . . . . . . . . 42 3.2.2 Physiological ranges of the physical model parameters 46 3.2.3 Local identifiability of the proposed state-space model 47 3.2.4 Grey-box identification . . . . . . . . . . . . . . . . . 50 3.2.5 Estimation of the bias and variance using a Monte- Carlo study . . . . . . . . . . . . . . . . . . . . . . . 54 3.2.6 Performance of the proposed grey-box identification method using simulated data . . . . . . . . . . . . . . 55 3.2.7 Estimation of the cardiac output from the radial artery pressure waveform . . . . . . . . . . . . . . . . . . . . 57 3.3 Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 vii Table of Contents 3.3.1 Data acquisition . . . . . . . . . . . . . . . . . . . . . 64 3.3.2 Data processing . . . . . . . . . . . . . . . . . . . . . 64 3.3.3 Cardiac output determination . . . . . . . . . . . . . 66 3.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 3.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 4 Measuring Aortic Pressure Waveform using Ultrasound B- Mode Imaging1,2 . . . . . . . . . . . . . . . . . . . . . . . . . . 75 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 4.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.2.1 Adaptive block matching . . . . . . . . . . . . . . . . 78 4.2.2 Gradient-based edge refining . . . . . . . . . . . . . . 80 4.2.3 Weighted total least-squares fitting . . . . . . . . . . 81 4.2.4 Continuity of the aortic wall motion . . . . . . . . . . 83 4.2.5 Choosing a cross-section along the arch . . . . . . . . 83 4.2.6 Estimating the pressure waveform from the distension waveform . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.3 Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 4.3.1 Data acquisition and processing . . . . . . . . . . . . 84 4.3.2 Calibration of the carotid artery pressure waveform obtained using applanation tonometry . . . . . . . . . 85 4.3.3 Comparison of waveforms . . . . . . . . . . . . . . . . 85 4.3.4 Calculations of pulse pressure, augmentation index and total arterial compliance . . . . . . . . . . . . . . 85 4.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.4.1 Waveform comparisons . . . . . . . . . . . . . . . . . 87 4.4.2 Comparisons of pulse pressures, augmentation indices and total arterial compliances . . . . . . . . . . . . . 87 4.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 4.6 Conclusion and future work . . . . . . . . . . . . . . . . . . . 94 5 Conclusions and Future Work . . . . . . . . . . . . . . . . . . 96 5.1 Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . 96 5.1.1 Model-based estimation of the aortic pulse wave ve- locity . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 5.1.2 Minimally-invasive monitoring of the cardiac output . 97 5.1.3 Estimation of the aortic pressure waveform using B- mode imaging . . . . . . . . . . . . . . . . . . . . . . 99 5.1.4 Suggestions for future work . . . . . . . . . . . . . . . 99 viii Table of Contents References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 Appendices A The Analogy between a Transmission Line and an Elastic Tube with Incompressible Flow . . . . . . . . . . . . . . . . . 121 B Implementation of the Transmission Line Model using the Finite Element Method . . . . . . . . . . . . . . . . . . . . . . 123 C The Derivation of the Transfer Function from the Ascending Aortic Flow to the Radial Artery Pressure . . . . . . . . . 125 ix List of Tables 2.1 Clinical characteristics of the subjects of the pulse wave ve- locity study . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.2 Fixed model parameters and initial model parameters . . . . 28 2.3 Identified model parameters using white box identification . . 29 2.4 Pulse wave velocity measured using different methods . . . . 29 2.5 Correlations between the echocardiographic pulse wave veloc- ity and different model-based pulse wave velocities . . . . . . 30 2.6 Performance of the transmission line model in predicting the aortic flow . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.1 Physiological ranges of the transmission line model parameters 48 3.2 The number of stable models for different time steps and num- bers of elements . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.3 The estimated bias and variances of the grey-box identifica- tion method . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.4 Model parameters estimated from the simulated data . . . . . 58 3.5 Stability of the models identified from the simulated data . . 59 3.6 Demographics of patients of the cardiac output monitoring project . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 3.7 Identified transmission line model parameters . . . . . . . . . 69 3.8 Cardiac output determinations . . . . . . . . . . . . . . . . . 70 3.9 Bias and limits of agreement in terms of % error . . . . . . . 70 3.10 Bias and limits of agreement in terms of cardiac output . . . 71 3.11 Bias and limits of agreement in terms of cardiac index . . . . 71 4.1 Demographics of subjects of the aortic pressure waveform project . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.2 Comparisons of the waveforms obtained using applanation tonometry and B-mode imaging in the frequency domain . . . 88 4.3 Comparisons of the derived pulse pressures, augmentation in- dices and total arterial compliances . . . . . . . . . . . . . . . 90 x List of Figures 1.1 The overall objectives . . . . . . . . . . . . . . . . . . . . . . 14 2.1 A single uniform transmission line terminated by a complex load . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 2.2 The calculation of the total arterial compliance . . . . . . . . 23 2.3 The model-predicted ascending aortic flow . . . . . . . . . . . 25 2.4 Noninvasive data collection of the model-based pulse wave velocity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 2.5 Box and whisker plots for the model-based pulse wave velocity and echocardiograhic pulse wave velocity . . . . . . . . . . . . 31 2.6 Correlation between the model-based pulse wave velocity and echocardiograhic pulse wave velocity . . . . . . . . . . . . . . 32 2.7 Bland-Altman plot of the model-based pulse wave velocity and echocardiograhic pulse wave velocity . . . . . . . . . . . . 33 2.8 A sample of the model-predicted ascending aortic flow wave- form and the corresponding experimentally-measured ascend- ing aortic flow waveform . . . . . . . . . . . . . . . . . . . . . 34 3.1 A lossless transmission line terminated by a complex load . . 40 3.2 The indexing for the current and the voltage in the finite element representation . . . . . . . . . . . . . . . . . . . . . . 43 3.3 The grey-box identification algorithm . . . . . . . . . . . . . 54 3.4 The simulated flow waveform used as the input to the trans- mission line model . . . . . . . . . . . . . . . . . . . . . . . . 60 3.5 The estimated radial artery pressure waveform using simu- lated flow as the model input . . . . . . . . . . . . . . . . . . 61 3.6 The re-estimation of the terminating resistance . . . . . . . . 61 3.7 The determination of the time constant of the Windkessel load using the area method . . . . . . . . . . . . . . . . . . . 62 3.8 The cardiac output estimation technique . . . . . . . . . . . . 65 3.9 The estimation of the left ventricular outflow tract flow from the Doppler spectrum . . . . . . . . . . . . . . . . . . . . . . 66 xi List of Figures 3.10 A sample of the radial artery pressure waveform . . . . . . . 67 3.11 A sample of the measured ascending aortic flow waveform . . 67 3.12 Log of singular value versus system order . . . . . . . . . . . 69 3.13 Covergence graphs . . . . . . . . . . . . . . . . . . . . . . . . 71 3.14 Bland-Altman plot of the cardiac output determinations . . . 72 3.15 Measured and estimated ascending aortic flow waveforms . . 73 4.1 Block selections for the block matching scheme . . . . . . . . 78 4.2 Search of the best-matching block . . . . . . . . . . . . . . . 79 4.3 The adaptive block matching scheme . . . . . . . . . . . . . . 80 4.4 Refining the edges using the gradient of the intensity . . . . . 81 4.5 Weight assigning scheme for fitting a second degree polyno- mial to the edge . . . . . . . . . . . . . . . . . . . . . . . . . 82 4.6 Comparisons of the waveforms obtained using applanation tonometry and B-mode imaging I . . . . . . . . . . . . . . . . 91 4.7 Comparisons of the waveforms obtained using applanation tonometry and B-mode imaging II . . . . . . . . . . . . . . . 92 xii Abbreviations and Symbols Chapter 2 AAF (cm3/s) . . . . . . . . . Ascending aortic flow AAP (mmHg) . . . . . . . . Ascending aortic pressure BSA (m2) . . . . . . . . . . . . Body surface area MFS . . . . . . . . . . . . . . . . . . Marfan syndrome PTT (ms) . . . . . . . . . . . . Pulse transit time PWV (cm/s) . . . . . . . . . Pulse wave velocity RMSE (cm3/s) . . . . . . . Root-mean-squared error TAC (cm3/mmHg) . . . Total arterial compliance DP (mmHg) . . . . . . . . . . Diastolic pressure at the brachial artery MP (mmHg) . . . . . . . . . . Mean arterial pressure SP (mmHg) . . . . . . . . . . Systolic pressure at the brachial artery C (cm3/mmHg/cm) . . Transmission line model unit line capacitance CL (cm 3/mmHg) . . . . . Transmission line model load capacitance CT (cm 3/mmHg) . . . . . Transmission line model total capacitance, sum of total line capacitance and load capac- itance L (mmHgs2/cm3/cm) . Transmission line model unit line inductance RL (mmHgs/cm 3) . . . . Transmission line model load resistance RT (mmHgs/cm 3) . . . . Transmission line model terminating resis- tance l (cm) . . . . . . . . . . . . . . . . Transmission line model line length Aao (cm 2) . . . . . . . . . . . . Ascending aortic cross-sectional area J (cm6/s) . . . . . . . . . . . . Performance index P (mmHg) . . . . . . . . . . . Pressure Q (cm3/s) . . . . . . . . . . . . Flow xiii Abbreviations and Symbols Qao,e (cm 3/s) . . . . . . . . . Experimentally-measured ascending aortic positive flow Qao,m (cm 3/s) . . . . . . . . Model-predicted ascending aortic positive flow ρb (mmHgs 2/cm2) . . . . Blood density µ . . . . . . . . . . . . . . . . . . . . . . Penalty coefficient of the performance index Λ (cm6/s) . . . . . . . . . . . . Penalty term of the performance index Chapter 3 AAF (cm3/s) . . . . . . . . . Ascending aortic flow AL (L) . . . . . . . . . . . . . . . Area integral of the model-estimated ascend- ing aortic flow waveform CO (L/min) . . . . . . . . . . Cardiac output D (min) . . . . . . . . . . . . . . Duration of the model-estimated ascending aortic flow waveform in one cardiac cycle ICU . . . . . . . . . . . . . . . . . . . Intensive care unit LVOT (cm) . . . . . . . . . . . Left ventricular outflow tract PAP (mmHg) . . . . . . . . . Peripheral artery pressure RAP (mmHg) . . . . . . . . Radial artery pressure SS. . . . . . . . . . . . . . . . . . . . . State-space TAC (cm3/mmHg) . . . Total arterial compliance C (cm3/mmHg/cm) . . Transmission line model unit line capacitance CL (cm 3/mmHg) . . . . . Transmission line model load capacitance Ctot (cm 3/mmHg) . . . . Transmission line model total line capaci- tance L (mmHgs2/cm3/cm) . Transmission line model unit line inductance Ltot (mmHgs 2/cm3) . . Transmission line model total line inductance RL (mmHgs/cm 3) . . . . Transmission line model load resistance RT (mmHgs/cm 3) . . . . Transmission line model terminating resis- tance Zch (mmHgs/cm 3) . . . . Transmission line model characteristic impedance of the line l (cm) . . . . . . . . . . . . . . . . Transmission line model line length Cwind (cm 3/mmHg) . . Windkessel model capacitance τ (s) . . . . . . . . . . . . . . . . . . Windkessel model time constant xiv Abbreviations and Symbols A(θ), B(θ), C(θ) . . . . . . Physical state-space model Â, B̂, Ĉ . . . . . . . . . . . . . . . Empirical state-space model T . . . . . . . . . . . . . . . . . . . . . Transformation matrix between the physical state-space model and the empirical state- space model e . . . . . . . . . . . . . . . . . . . . . . Error between the physical state-space model and the empirical state-space model A1, A2 . . . . . . . . . . . . . . . . Parameters of the state-space transmission line model Ij . . . . . . . . . . . . . . . . . . . . . Current at the jth node of the finite element transmission line K1, K2, K3, K4 . . . . . . . Parameters of the state-space transmission line model M . . . . . . . . . . . . . . . . . . . . . Transfer function from the ascending aortic flow to the radial artery pressure Vj . . . . . . . . . . . . . . . . . . . . . Voltage at the jth node of the finite element transmission line d . . . . . . . . . . . . . . . . . . . . . . Discrete-time delay from the ascending aortic flow to the radial artery pressure n . . . . . . . . . . . . . . . . . . . . . . Discrete time variable p− 1 . . . . . . . . . . . . . . . . . . Number of elements of the finite element transmission line t . . . . . . . . . . . . . . . . . . . . . . Continuous time variable x . . . . . . . . . . . . . . . . . . . . . . Continuous spatial variable Φc . . . . . . . . . . . . . . . . . . . . Voltage across CL ∆t . . . . . . . . . . . . . . . . . . . . Time step α . . . . . . . . . . . . . . . . . . . . . . Parameter of the state-space transmission line model β . . . . . . . . . . . . . . . . . . . . . . Correction factor between the model-estimated cardiac output and echocardiograhic-measured cardiac out- put ϕ . . . . . . . . . . . . . . . . . . . . . . Vector of unknown physical model parame- ters µ . . . . . . . . . . . . . . . . . . . . . . Propagation velocity of the transmission line θ . . . . . . . . . . . . . . . . . . . . . . Vector of unknown parameters of the state- space transmission line model ψ . . . . . . . . . . . . . . . . . . . . . . Vector of known parameters of the state- space transmission line model xv Abbreviations and Symbols Chapter 4 AI (%) . . . . . . . . . . . . . . . Augmentation index PP (mmHg) . . . . . . . . . . Pulse pressure RMSD (mmHg) . . . . . . Root-mean-squared difference TAC (cm3/mmHg) . . . Total arterial compliance DP (mmHg) . . . . . . . . . . Diastolic pressure at the brachial artery MP (mmHg) . . . . . . . . . . Mean arterial pressure SP (mmHg) . . . . . . . . . . Systolic pressure at the brachial artery A(t) (cm2) . . . . . . . . . . . . Time function of the cross-sectional area of the ascending aorta A0 (cm 2) . . . . . . . . . . . . . Cross-sectional area at the reference pressure D (cm) . . . . . . . . . . . . . . . Diameter of the ascending aorta P (t) (mmHg) . . . . . . . . . Time function of the ascending aortic pres- sure P0 (mmHg) . . . . . . . . . . . Reference pressure Pa (mmHg) . . . . . . . . . . . Pressure calculated from the exponential pressure-area relationship Pdw (mmHg) . . . . . . . . . Pressure waveform obtained from the disten- sion waveform Pexp (mmHg) . . . . . . . . Pressure waveform obtained from the disten- sion waveform using exponential calibration Pinf (mmHg) . . . . . . . . . Pressure at the inflection point Plin (mmHg) . . . . . . . . . Pressure waveform obtained from the disten- sion waveform using linear calibration Pmax (mmHg) . . . . . . . Maximum pressure Ptono (mmHg) . . . . . . . Pressure waveform obtained using applana- tion tonometry RP (mmHgs/cm 3) . . . . Peripheral resistance cr . . . . . . . . . . . . . . . . . . . . . Correlation between the kernel block and the candidate block in block matching cv . . . . . . . . . . . . . . . . . . . . . Curvature of the second degree polynomial fit k . . . . . . . . . . . . . . . . . . . . . . Size of the kernel block in block matching s . . . . . . . . . . . . . . . . . . . . . . Size of the search region in block matching vx . . . . . . . . . . . . . . . . . . . . . x-coordinate of the second degree polynomial fit xvi Abbreviations and Symbols vy . . . . . . . . . . . . . . . . . . . . . y-coordinate of the second degree polynomial fit β . . . . . . . . . . . . . . . . . . . . . . Threshold of discrepancy in the gradient- based edge refining scheme ϕ . . . . . . . . . . . . . . . . . . . . . . Threshold of distance in the weighted total least-squares fitting scheme γ . . . . . . . . . . . . . . . . . . . . . . Threshold of correlation in the adaptive block matching scheme θ . . . . . . . . . . . . . . . . . . . . . . Threshold of the number of neighborhood points in the weighted total least-squares fit- ting scheme xvii To my mum xviii Chapter 1 Introduction 1.1 Motivation According to the World Health Organization (WHO), cardiovascular dis- eases are the number one cause of death globally, and are projected to re- main the single leading cause of death by 2030. Central pressure and arterial stiffness have strong association with the pathogenesis of the cardiovascular diseases [1]. Central arteries are proximal to the heart and include the aorta and the carotid artery while peripheral arteries are distal to the heart and include the arteries at the limbs, such as the radial, brachial and femoral arteries. Although brachial pulse pressure has long been recognized as a powerful predictor of cardiovascular events [103,113], it has recently been shown that central pulse pressure has stronger association with the severity of coronary artery disease than brachial pulse pressure [168], and is a powerful indepen- dent predictor of all-cause, including cardiovascular, mortality in end-stage renal diseases [136]. In addition, it has also been suggested that central pres- sure may have predictive value independent of the corresponding brachial pressure [1]. Therefore, there is an increasing interest in the measurement of central pressure. However, central hemodynamics are more difficult to assess than the peripheral hemodynamics. Aortic pressure can be measured by advancing a catheter through the brachial artery or the femoral artery to the aorta, but such procedure is highly invasive. In contrast, pressure at the superficial arteries, including carotid, radial and femoral arteries, can be measured non-invasively using applanation tonometry. It has been suggested that carotid artery pressure obtained using appla- nation tonometry can be used as a surrogate of aortic pressure [82]. However, the procedure requires substantial expertise. In addition, the aortic pressure and the carotid artery pressure are not identical due to difference in location, and the difference may become significant in case of carotid artery diseases, such as carotid stenosis. In this thesis, a new method is developed to mea- sure the aortic pressure using B-mode imaging. This method estimates the 1 1.1. Motivation aortic pressure waveform using B-mode images of the aortic arch, which can be obtained during routine echocardiographic exam, and therefore does not require any special equipment or training. It is applicable to patients with carotid diseases. Since this method is non-invasive, it can be used in early detection of cardiovascular diseases in the general population. Aortic pulse wave velocity is the most well-accepted index of arterial stiffness. It has been shown that aortic pulse wave velocity is an independent predictor of all-cause and cardiovascular mortality in hypertensive patients [92], cardiovascular mortality in subjects over 70 years of age [109], and cardiovascular events in different populations [19,93,154]. The most common method to estimate the aortic pulse wave velocity is to measure the pressure pulses at the carotid and femoral arteries, and calculate the aortic pulse wave velocity as the distance between the measuring sites divided by the pulse transit time. This method requires high precision of the measurement of the pulse transit time. The distance is measured on the surface of the skin and may be inaccurate in subjects with large busts or abdominal obesity. In this thesis, a non-invasive, model-based method is developed to measure the aortic pulse wave velocity. This new method does not require pulse transit time or distance measurement. Since this method is non-invasive, it can be used in early detection of cardiovascular diseases in the general population. The modelling approach has been used by previous researchers to esti- mate the stiffness of the central arteries. The most commonly used model to represent the arterial system is the two-element Windkessel [119]. The two-element Windkessel model is simple and captures the two most impor- tant mechanical properties of the arterial system: compliance and resistance. However, the Windkessel model is a lumped model and has two related lim- itations. First, it implies that pulse wave velocity is infinite. Second, it cannot capture the wave propagation and reflection characteristics in the arterial system. Wave reflection is an important phenomenon in the arterial system and is related to the pathogenesis of the cardiovascular diseases [1]. The transmission line model is distributive and is a more suitable model for the arterial system [94]. It is one of the motivations of this thesis to use the more accurate transmission line model to estimate the stiffness of the central arteries. Non-invasive methods to measure the central pressure and arterial stiff- ness are important as they are valuable tools for the diagnosis of cardiovas- cular diseases at their early stage of development. Invasive measurement methods are associated with significant risk and are not always justifiable. Non-invasive methods are also useful for large-scale, cross-sectional studies. Cardiac output indicates the performance of the heart as a pulsatile 2 1.2. Background pump [127]. It has been shown that the cardiac index in postoperative patients is associated with their survival rate [18,40]. The monitoring of cardiac output is therefore essential during perioperative care. Cardiac out- put can be measured using thermodilution, which is often performed with a pulmonary artery catheter. The placement of a pulmonary artery catheter is associated with increased mortality and increased utilization of resources [34]. A less invasive alternative is much desired. Cardiac output can be measured using transthoracic Doppler ultrasound, but the procedure requires a well-trained operator, and therefore is not suitable for long-term continuous monitoring. In contrast, the radial artery pressure can be monitored continuously using an indwelling catheter and the procedure is minimally invasive. Therefore, in this thesis, a new technique is developed to monitor the cardiac output using the radial artery pressure waveform. This technique can be applied to post-operative patients during their stay at the intensive care unit. The estimation of the cardiac output using the radial artery pressure waveform is termed pulse contour analysis. One of the most popularly used pulse contour analyses is Wesseling’s algorithm [172], which is based on the a simplified transmission line model. According to Wesseling’s transmission line model, the aortic flow and the aortic pressure are directly proportional and are related by the characteristic impedance of the line. However, aortic flow is instead directly proportional to the difference between the forward aortic pressure and the backward aortic pressure. As presented in Section 1.3.3, Wesseling’s algorithm uses only the area of the systolic portion of the pressure waveform to estimate the cardiac output. In addition, in Wessel- ing’s algorithm, the radial artery pressure is used as a surrogate of aortic pressure. It is one of the motivations of this thesis to develop a cardiac out- put estimation technique based on a more precise transmission line model. 1.2 Background 1.2.1 The mechanism of the pressure and flow waves in the arterial system The cardiovascular system can be divided into two major components, the heart and the vasculature. While the heart is a muscular pump, the vascu- lature is a network of vessels, which consists of the arteries, the capillaries and the veins. Blood runs from the ventricles of the heart, through the ar- teries, the capillaries and the veins, and back to the atria of the heart. Since the heart operates as an intermittent pump, the pressure and flow waves 3 1.2. Background entering the arteries have a oscillating component and a steady component. The oscillating component is gradually eliminated in the arterial system, and only the steady component enters the capillaries. The arteries form a network of progressively branching, geometrically and elastically tapering viscoelastic tubes. The arterial network has two major functions, cushioning and conduiting. The arteries are compliant and are able to distend to accommodate the sudden increase in blood volume, thereby providing a cushioning function. The arteries then recoil slowly and push the volume of blood forward. There is a progressive loss of the cushioning function along the arterial network. The aorta and its major branches are referred to as the elastic arteries. The media of the elastic arteries are rich in elastic fibers. The elastic arteries gradually merge into the muscular arteries, and the muscular arteries gradually merge into the arterioles. The media of the muscular arteries is dominated by smooth muscle cells and has very few elastic fibers. From the elastic arteries to the resistive arterioles, the compliance of the arteries decreases gradually [119]. Pressure and flow waves propagate down the arterial network at a finite pulse wave velocity. During the propagation, wave reflections occur at var- ious discontinuities and terminations, including the branching points and the resistive arterioles [121]. Theoretically, since the compliance and the cross-sectional area of the arteries are progressively changing, wave reflec- tions occur at every point along the arterial network. However, it is believed that the major wave reflections occur at the resistance arteries, including the muscular arteries and the arterioles [1]. A direct consequence of wave re- flections is pressure amplification from the central arteries to the peripheral arteries. The pressure wave recorded at any point along the arterial network is the sum of the incident wave and all the reflected waves. At the periph- eral sites, the reflected wave augments the systolic portion of the incident wave. The reflected wave propagates from the peripheral site back to the central site and is gradually attenuated. At the central sites, the attenuated reflected wave augments the diastolic portion of the incident wave. The net result is that the pressure amplitude is higher at the peripheral site than the central site. This phenomenon is referred to as pressure amplification. During systole, it is the aortic systolic pressure that the left ventricle encounters. Chronic sustained elevation of aortic systolic pressure causes left ventricular hypertrophy [122]. In young and healthy subjects, reflected wave augments the incident wave at the ascending aorta during diastole. In the case of arteriosclerosis, which is the degeneration of the arterial media, the arterial stiffness increases. According to the Moen- Korteweg equation, pulse wave velocity increases as arterial stiffness increases [119]. Because of 4 1.3. Literature review the increase in pulse wave velocity, the reflected wave arrives at the ascend- ing aorta during systole and causes the aortic systolic pressure to increase. Arteriosclerosis occurs mainly in aging and hypertensive subjects [122]. Sim- ilarly, in the case of arterial or arteriolar constriction, where the reflection site is moved closer to the heart, the reflected wave augments the incident wave at the ascending aorta during systole causing aortic systolic pressure to increase. In summary, the central pulse pressure, the aortic pulse wave velocity and the central arterial stiffness are interrelated and have strong influences on the pathogenesis of the cardiovascular diseases. Because of pressure amplifi- cation, peripheral pressure cannot be used a as surrogate of central pressure. 1.3 Literature review 1.3.1 Current non-invasive methods to measure the central pressure waveform Carotid artery applanation tonometry Arterial applanation tonometry is a technique for transcutaneous pressure recording. It is based on the principle that the force required to applanate an artery divided by the applanated area is equivalent to the intralumenal pressure. The principle of applanation tonometry was originally applied in ocular tonometry to measure the intraocular pressure. In 1963, Pressman and Newgard published the earliest work in arterial applanation tonometry [129]. Unlike sphygmomanometry, arterial applanation tonometry can mea- sure the instantaneous blood pressure continuously and thereby provide the pressure waveform. The early development of arterial applanation tonome- try was described in the review by Drzewiecki et al. [43]. In 1989, Kelly et al. [80] developed a tonometer consisting of a pencil-type probe incorporating a high-fidelity Millar gauge transducer, and it was shown that the tono- metrically and intra-arterially measured radial artery pressure waveforms were similar. It was also shown that the morphology of the tonometrically measured carotid artery pressure waveform and that of the intra-arterially measured ascending aortic pressure waveform were similar [81]. During the early development of arterial applanation tonometry, re- searchers aimed to measure the absolute pressure. It is then regarded as possible only if stringent requirements of applanation are met [81]. External calibration can be applied to the brachial artery pressure waveform. The brachial systolic and diastolic pressures are measured at the arm using cuff 5 1.3. Literature review sphygmomanometer. The systolic pressure is assigned to the peak and the diastolic pressure is assigned to the trough of the waveform. For the ra- dial artery, it is often assumed that the systolic and diastolic pressures at the brachial artery and the radial artery are the same. The same calibra- tion is thus applied. However, because of pressure amplification, the same calibration scheme cannot be applied to the carotid artery. Kelly et al. [82] proposed an external calibration method for the carotid artery. It is based on the assumption that the mean arterial pressure and the diastolic pressure remain constant through out the arterial tree. The brachial artery pressure waveform is measured using applanation tonometry and calibrated using a cuff sphygmomanometer. The mean arterial pressure is calculated as the area integral over the duration of the calibrated brachial artery pressure waveform. The mean arterial pressure is then assigned to the area integral over the duration of the carotid artery pressure waveform. The diastolic pressure is assigned to the trough of the carotid artery pres- sure waveform. It was shown that the calibrated, tonometrically measured carotid artery pressure waveform and the inter-arterially measured aortic pressure waveform were similar [82]. Carotid artery applanation tonometry has since become a popular choice to measure the central pressure waveform [7, 65, 101, 138, 167]. In some studies, a simplified calibration process was used [3, 128, 134]. The mean arterial pressure (MP) was calculated from the systolic pressure (SP) and diastolic pressure (DP) at the brachial artery using the traditional formula: MP = 13SP + 2 3DP instead of the brachial artery pressure waveform. The major limitation of carotid artery applanation tonometry is that the applanation of the carotid artery is difficult and requires substantial amount of expertise. The quality of the measured waveform is affected by the lateral and radial motions of the applanation force. The carotid artery is surrounded by loose tissue and is therefore not an ideal site for optimal applanation. Radial artery applanation tonometry with a generalized transfer function Comparing to the carotid artery, the radial artery is well-supported by bony tissue allowing optimal applanation to be achieved. Because of pressure amplification, the pressure waveforms at the ra- dial artery and the central arteries are very different in both morphology and magnitude. Radial artery pressure waveform can not be used as a surrogate of the central pressure waveform. Karamanoglu et al. [77] pro- 6 1.3. Literature review posed to use a generalized transfer function to reconstruct the aortic pres- sure waveform from the radial artery pressure waveform. The generalized transfer function is calculated as the average of the individual transfer func- tions from the ascending aortic pressure to the radial artery pressure of a group of subjects. Since then, the idea was explored by a number of studies [32,46,68,69,126]. Most of these studies derived a non-parametric transfer function [32,68,69] while one study also derived an ARX transfer function [46]. In these studies, the radial artery pressure waveform was measured either intra-arterially [126], or tonometrically with invasive calibration, in which the the mean and diastolic pressures of the radial artery and the aorta were matched [32,46,68,69]. This is to eliminate the error introduced by the cuff sphygmomanometer during the calibration of the radial artery pressure waveform. Three of these studies concluded that the generalized transfer function method is viable [32,46,126], while two were skeptical [68,69]. In the study by Hope et al. [68], it was shown that even though the mean er- ror between the invasively-measured and the transfer-function-derived aortic systolic pressures was small, the error variance was large. It was pointed out that a generalized transfer function derived from a population mean would be expected to result in a small mean error in the estimation of the systolic pressure or other parameters, however, this does not demonstrate its validity on an individual basis. The large error variance indicated that the generalized transfer function may not be valid for each individual. A commercial device, SphygmoCor (SphygmoCor; Atcor Medical, Syd- ney, Australia), estimates the aortic pressure waveform using the generalized transfer function method. A substantial amount of studies used this device to measure central pulse pressure and other central indices. However, the validity of the generalized transfer function method remains controversial [70]. A number of studies have shown that SphygmoCor cannot measure central pulse pressures with acceptable accuracy [33,37,147,149]. Distension waveform method The relationship between the pressure waveform and the distension wave- form of an artery is determined by the elastic properties of the arterial wall. In 1978, Ghista et al. [52] published the earliest theoretical analysis for the determination of the pressure waveform from the distension waveform. It was suggested that if several elasticity parameters are known, the pressure waveform can be calculated from the distension waveform. It was also in- dicated that high resolution for the measurement of the arterial diameter is required. The idea was not explored further for about twenty years, and it 7 1.3. Literature review was most likely due to the complexity in the determination of the elasticity parameters and the lack of accuracy in the measurement of the distension waveform. In 1998, Watanabe et al. [171] used an echo-tracking device to measure the aortic distension waveform, and a relatively simple exponential relationship between the pressure and the diameter to calculated the aortic pressure waveform. The echo-tracking device can provide a spatial resolu- tion up to 20 µm. However, the simple exponential relationship involves a stiffness parameter, and the study did not provide a method to accurately estimate this parameter. In 2001, Van Bortel et al. [161] used the echo- tracking device to measure the carotid distension waveform, and assumed a linear relationship between the pressure and the diameter. With carotid artery applanation tonometry as the reference method, the central pulse pressure obtained from the carotid distension waveform demonstrated good accuracy. In 2004, Meinders and Hoeks [110] used high frame rate B-mode images to measure the carotid distension waveform, and an exponential re- lationship between the pressure and the cross-section area of the artery to calculate the pressure waveform . The exponential relationship involves a stiffness parameter. An iterative scheme was proposed to find the stiffness parameter. It was assumed that the mean pressures at the brachial artery and the carotid artery were the same. The stiffness parameter was updated iteratively such that the mean pressure of the calibrated distension wave- form matched that of the brachial artery. In 2008, Vermeersch et al. [166] used the echo-tracking device to measure the carotid distension waveform, and both linear and exponential calibration schemes to calculate the pres- sure waveform. It was shown that the exponential calibration scheme was superior to the linear calibration scheme, and good accuracy was archived in the estimation of the carotid artery pressure waveform and pulse pressure. Compared to applanation tonometry, the distension waveform method is relatively new and is not yet popularly used in clinical research for the measurement of the central pressure waveform. 1.3.2 Current non-invasive methods to measure aortic pulse wave velocity Pulse wave velocity is often measured using the two-site method. The pulse wave is measured at two recording sites, and the pulse wave velocity is cal- culated as the distance between the sites divided by the time delay between the measured pulse waves. The time delay is often referred to as the pulse transit time. Pulse wave velocity can be obtained from the travel of the pres- sure wave, the flow wave or the distension wave. Aortic pulse wave velocity 8 1.3. Literature review is most commonly measured as the carotid-femoral pulse wave velocity. In 1922, Bramwell and Hill [23] published the earliest clinical study on the measurement of the pulse wave velocity. The pressure waveforms at the carotid artery and the radial artery were measured simultaneously us- ing the hot-wire sphygmograph. The pulse transit time was taken as the time interval elapsed between the commencements of the carotid and ra- dial deflections displayed on the sphygmogram. The distance traveled was taken as the distance from the sterno-clavicular joint to the wrist, less the distance from the same point to the carotid artery. It should be noted that the pressure waveforms obtained using the sphygmograph in this study were quite different from the intra-arterially measured pressure waveforms in the literature. However, large deflections corresponded to the systolic upstrokes can be easily identified, and the pulse transit time was measured from the deflections. It should also be noted that this method measures the pulse wave velocity along the subclavian-radial pathway and does not measure the aortic pulse wave velocity. In contrast to the carotid-radial pulse wave velocity, the carotid-femoral pulse wave velocity measures the pulse wave velocity along the aorta-iliac pathway and is a direct measurement of the aortic pulse wave velocity. Several studies measured the carotid-femoral pulse wave velocity. Some of these studies used the flow waves obtained using Doppler ultrasound [5,8,9,100,101]. Some of these studies used the pulse waves obtained using a pulse transducer [4,73,124,135]. It was not clarified if the pulse transducer measured the pressure waves or the flow waves. However, it is likely that the pulse wave measured was related to the pressure wave as the studies also mentioned a diastolic curve as part of the pulse wave [4,124,135]. A study used the pressure wave obtained using a pressure probe [90]. The principle of the pressure probe was not explained nor was the quality of the measured pressure wave displayed. A study used the distension wave obtained by echo-tracking [162]. The study used an echo-tracking device evaluated by Kool et al [88], and it was shown that distension waves with close resemblance of the typical pressure waves can be obtained using this device. Currently, two commercial devices, SphygmoCor and Complior, are pop- ularly used to obtain the carotid-femoral pulse wave velocity in clinical re- search. Using the SphygmoCor device, the pressure waveforms at the carotid artery and the femoral arteries are acquired consecutively using applanation tonometry [45,83,111]. The pressure waveforms are synchronized using the electrocardiographic signal as the reference. The pulse transit time is cal- culated as the foot-to-foot interval between the pressure waveforms. The 9 1.3. Literature review distance between the sites is taken as the distance from the sternal notch to the applanation site at the femoral artery as measured on the surface of the skin [83,111], or the distance from the sternal notch to the femoral applanation site less the distance from the sternal notch to the carotid ap- planation site [45]. Using the Complior device, the pressure waveforms at the carotid artery and femoral arteries are acquired simultaneously using a piezo-electric pressure transducer [6,10,17,104]. The pulse transit time is determined by a correlation algorithm. In short, the correlation between the part of the initial rise to just after the peak of the carotid waveform and that of the time-shifted femoral waveform is found, and the pulse transit time is estimated from the time-shift that gives the maximum correlation [6]. The distance is taken as the distance between the two applanation sites [6,10,104]. In the study by Asmar et al [6], the pressure waveforms obtained using the pressure transducer were shown, and they were not similar to the typical pressure waveforms found in literature. However, clear initial upstrokes of the waveforms can be easily identified. SphygmoCor and Complior use different algorithms to measure the pulse transit time and the distance traveled. It has been shown that the pulse transit time calculated using these two algorithms are significantly different and the measurement of pulse wave velocity is critically dependent on the algorithm used to measure the pulse transit time [112]. It has also been shown that travel distances measured by these two methods are significantly different [130]. In addition, abdominal obesity and large busts in women can also make distance measurement inaccurate. Pulse wave velocity can also be measured along the aorta using MRI [58,133,141]. Flow velocity is obtained at different levels along the aorta using velocity-encoded MRI images [58]. The pulse transit time is obtained from the velocity curves. Distances between the levels are measured on the long-axis image of the aorta [58]. Pulse wave velocity can also be measured along the aortic arch using Doppler ultrasound [22, 140]. Flow waveforms are measured at the ascending and descending aorta. The pulse transit time is calculated as the foot-to-foot interval between the measured flow waveforms. The distance is measured on the B-mode image of the aortic arch. Recently, a new device, the arteriograph, has been developed to measure aortic pulse wave velocity [11,75]. This device makes use of an oscillatory method to measure the brachial artery pressure waveform through an in- flated cuff. In short, the change in the pressure in the artery beneath the inflated cuff induces oscillations, and these oscillations are detected and used as an indirect measure of the change in the pressure. To measure the pulse 10 1.3. Literature review wave velocity, the cuff is over-inflated to 35 mmHg over the systolic pressure. It is claimed that, because of the over-inflation, the signal of the reflected wave is intensified and can be detected more easily [75]. The time interval between the first wave and the second wave of the measured waveform cor- responds to the distance from jugulum to the symphysis, and the pulse wave velocity is thus derived [11]. It has been shown that the pulse wave velocity measured by this device is related to the pulse wave velocities measured by SphygmoCor and Complior [11,75] and also invasively [71].However, it has also been shown that the relation between the pulse wave velocity measured by arteriograph and the carotid-femoral pulse wave velocity measured by applanation tonometry is low [118]. In addition, the underlying principle of arteriograph is questioned. Through a simulation study, it has been shown that arteriograph can only pick up the wave reflection phenomenon confined to the brachial artery and reflect the stiffness of the brachial artery [159]. 1.3.3 Current methods to monitor cardiac output Pulse contour analysis Pulse contour analysis is a method for the estimation of the cardiac output using the arterial pressure waveform. The idea was put forward as early as 1904 [62]. The central aorta is considered as a elastic chamber, which dis- tends during systole to accommodate the sudden increase in volume. There- fore, the rise in pressure during systole is related to the ventricular ejection. This relationship is dependent on various factors, such as the size and elas- ticity of the aorta, and therefore varies among individuals. Calibration is thus required to determine the absolute value of the cardiac output from the arterial pressure. Various equations have been put forward to describe the pressure-flow relationship [2,170]. Wesseling et al. [172] proposed to relate the area under the systolic portion of the pressure waveform to the cardiac output through the charac- teristic impedance of the aorta. The characteristic impedance is determined during calibration, and is updated during the estimation of the cardiac out- put according to various factors, such as the changes in the heart rate and the mean arterial pressure. Several studies have shown that this method can estimate cardiac output accurately with thermodilution as the refer- ence method [53,74,132,158]. However, it has also been suggested that this method may become unreliable during hemodynamic instability [132,158]. It should also be noted that, theoretically, aortic pressure waveform should be used as the input of Wesseling’s algorithm [74], but in practice the radial 11 1.3. Literature review or femoral artery pressure waveform is used as a surrogate. Currently, three commercial devices, PiCCO, PulseCO and FloTrac, es- timate cardiac output using pulse contour analysis. PiCCO is based on Wesseling’s algorithm [53,132]. The more recent ver- sion also uses the diastolic portion of the pressure waveform to estimate arterial compliance [48]. It has been shown that the enhanced version of PiCCO can estimate cardiac output with acceptable accuracy during hemo- dynamic instability [49]. PulseCO is based on a cardiac output-radial artery pressure relationship developed in the frequency domain [97]. To use this relationship to derive the cardiac output, three quantities, including the first harmonic of the aortic input impedance, the first harmonic of the transfer function from the aortic pressure to the radial artery pressure, and the ratio of the cardiac output to the first harmonic of aortic flow, have to be approximated. These three quantities are estimated from the radial artery pressure using equations that are derived from the data of the modelling studies by Karamangoglu et al. [78,79] and the in-vitro study by Langewouters et al. [91]. It has been shown that the limits of agreement between PulseCO and lithium dilution during a 4-hour calibration-free period was about 20% to 30%, which is acceptable [29]. Low bias but wide confidence limits of PulseCO has also been report with thermodilution as the reference method [108]. It has also been shown that PulseCO cannot estimate cardiac output reliably during change of systemic vascular resistance [176]. In contrast to the PiCCO system and the PulseCO system, the Flo- Trac/Vigileo system does not require external calibration. Instead, cardiac output is estimated from the arterial pressure waveform in conjunction with the demographic data of the patients. The relationship of the cardiac out- put and the arterial pressure waveform is based on the data from an in-vitro study of the elasticity of human aorta [105]. The validation of the accuracy of this device has been disappointing so far. In an article by Biancoflore et al. [16], the validation studies of this device published from 2006 to 2008 were reviewed. Only two out of 13 published studied achieved an agreement error of less than 30 %. Recently, it is claimed that a software upgrade can improve the performance of this device [106, 146]. However, studies have shown that even with the upgraded softwere, the agreement error of this device is still larger than 30 % [16,146]. Further investigation is required. 12 1.3. Literature review Transcutaneous echocardiography A commercial device, USCOM, estimates cardiac output using transcuta- neous echocardiography [157]. The device consists of a Doppler flow probe. The probe is placed on the chest at either the left parasternal position to measure pulmonary blood flow velocity or the suprasternal position to mea- sure aortic blood flow velocity. Stroke volume is calculated as the velocity- time integral of the measured Doppler spectrum times the cross-sectional area of the valve. The cross-sectional area is calculated from the subject’s height. Cardiac output is calculated as the product of the stroke volume and the heart rate. Training is required to use this device, but it has been suggested that the learning curve is short [35,157,164]. It has also been sug- gested that inter-observer variability may be high [125,157,]. Some studies have shown that the device has clinically acceptable agreement with ther- modilution [31,35] while some studies have shown poor agreement [20,30,87]. It has also been shown that the device has poor agreement with conventional echocardiography method [125]. It has also been suggested that more accu- rate measurement of the cross-sectional area may improve the performance of the device [87]. Other methods Cardiac output can be estimated using transoesophageal echocardiography. The instrument consists of a standard echocardiographic system attached to a conventional flexible endoscope. The distal end of the endoscope is advanced into the esophagus and the ultrasound transducer is positioned behind the heart. Sedation is usually provided to reduce the subject’s anx- iety and improve the subject’s compliance [50]. The procedure is usually considered as semi-invasive. Cardiac output can be estimated non-invasively using the partial carbon dioxide rebreathing method based on the Fick’s principle. In the direct Fick method, the mix venous oxygen or carbon dioxide content is measured from the blood sample taken at the pulmonary artery. In the rebreathing method, the subject breathes into a bag of gas mixture for a couple of times. The mixed venous carbon dioxide oxygen is then estimated from carbon dioxide concentration of the bag. The detailed description of the procedure is given by Franciosa [51]. This method is only accurate during precisely controlled mechanical ventilation and showed poor agreement with thermodilution during spontaneous breathing [156]. 13 1.4. Thesis objectives Central hemodynamics Pressure Flow Pulse wave velocity Non-invasive estimation of the aortic pressure waveform using B-mode ultrasound Minimally-invasive monitoring of the cardiac output Non-invasive measurement of the aortic pulse wave velocity Figure 1.1: The overall objectives 1.4 Thesis objectives The overall objective of this thesis is to develop and evaluate non-invasive and minimally-invasive methods for measuring and monitoring three aspects of central hemodynamics, including pressure, flow and the propagation ve- locity of pressure and flow pulses, which is equivalent to PWV (Figure 1.1). In this thesis, novel methods are developed to • measure the arterial stiffness using non-invasive data. • monitor the cardiac output using minimally-invasive data. • measure the aortic pressure waveform non-invasively. From a modelling perspective, the underlying algorithms of these methods are interrelated. The arterial system is modeled as a single uniform lossless transmission line terminated by a complex load. Each measuring or mon- itoring task is a problem of identification, simulation or data acquisition. The specific objectives of this thesis are as follows: 1. Formulating the arterial stiffness estimation problem as an identifica- tion problem. (a) Modelling the arterial system as a single lossless uniform trans- mission line. 14 1.4. Thesis objectives (b) Identifying the model using the carotid artery pressure and the aortic flow as the input and output to the model. (c) Estimating the arterial stiffness using the identified model pa- rameters. 2. Formulating the cardiac output monitoring problem as a problem of identification and simulation. (a) Modelling the pulse transmission path from the ascending aorta to the radial artery as a single lossless uniform transmission line. (b) Identifying the model using the aortic flow and the radial artery pressure as the input and output to the model. (c) Deriving a transfer function from the radial artery pressure to the aortic flow using the identified model parameters. (d) Simulating the aortic flow using the derived transfer function and the radial artery pressure. 3. Measuring the aortic pressure waveform using the aortic distension waveform obtained using B-mode imaging. (a) Developing a semi-automatic algorithm for the measurement of the aortic distension waveform from a series of B-mode images of the aortic arch. (b) Applying a calibration scheme to estimate the aortic pressure waveform from the aortic distension waveform. The following contributions were made: 1. As presented in Chapter 2, during the development of the arterial stiffness estimation method: (a) A white-box identification algorithm was developed for the trans- mission line model implemented using a finite-element method. (b) The white-box identification algorithm was applied to the non- invasive data obtained from a group of healthy children and a group of children with Marfan syndrome. Aortic pulse wave ve- locity and the ratio of the large artery compliance to the small artery compliance were obtained from the identified model pa- rameters. These two indices could differentiate the group of healthy children from the group of children with Marfan syn- drome. 15 1.4. Thesis objectives 2. As presented in Chapter 3, during the development of the cardiac output monitoring method: (a) A state-space representation of the transmission line model based on the finite element method was developed. (b) A grey-box identification algorithm was developed for the pro- posed state-space representation of the transmission line model. (c) The model identifiability of the proposed state-space representa- tion, and the bias and variances of the proposed grey-box identi- fication algorithm were investigated. (d) A cardiac output monitoring technique was developed using the proposed grey-box identification. The cardiac output monitor- ing technique was applied to a group of post-operative infants. The technique achieved clinically acceptable accuracy without re-calibration for the first five hours. 3. As presented in Chapter 4, during the development of the method for aortic pressure waveform measurement: (a) A semi-automatic algorithm for the measurement of the aortic distension waveform based on block-matching, edge detection us- ing the gradient of the intensity, and smoothness constraint on the continuity of the wall motion was developed. (b) The algorithm was applied to the B-mode images obtained from a group of healthy children. The aortic pressure waveform was estimated from the aortic distension waveform through linear and exponential calibration schemes, and compared to the carotid artery pressure waveform obtained using applanation tonometry. The morphologies of the waveforms, and the derived central pulse pressure and total arterial compliance were in good agreement. 16 Chapter 2 Model-based Aortic Pulse Wave Velocity1,2,3 Aortic pulse wave velocity (PWV) is an independent predictor of all-cause and cardiovascular mortality in hypertensive patients [92], cardiovascular mortality in subjects over 70 years of age [109], and cardiovascular events in different populations [19,93,154]. Aortic PWV is often estimated as the carotid-femoral pulse wave veloc- ity, which is measured as the pulse transit time (PTT) from the carotid site to the femoral site divided by the distance between the two sites [112,130]. This method is very sensitive to PTT. Different methods to measure PTT can result in significantly different PWV [112]. In addition, the distance be- tween the sites is measured on the surface of the skin, and may be inaccurate in subjects with large busts or abdominal obesity. Aortic PWV is also popularly measured using echocardiography [22,140]. PTT is determined from the Doppler flow pulses measured at the ascending and descending aorta, and the distance is measured on the image of the aortic arch. However, in some subjects, a clear image of the aortic arch may be difficult to obtain. To avoid the difficulties caused by PTT and distance measurements, a model-based method for estimating PWV without PTT and distance mea- surements is developed in this chapter. 1 A version of this chapter has been submitted. 2 M. Leung, G. A. Dumont, G. G. S. Sandor and J. E. Potts, “Estimating arterial stiffness using transmission line model,” 28th IEEE Engineering in Medicine and Biology Conference, New York, 2006 3 M. Leung, G. A. Dumont, J. E. Potts, M. T. Potts and G. G. S. Sandor, “A model- based method to estimate pulse wave velocity,” 11th World Congress of Medical Physics and Biomedical Engineering, Munich, Germany, 2009 17 2.1. Introduction 2.1 Introduction Arterial stiffness affects the change in the lumen area of the arteries given a change in pressure. It has been suggested that increase in arterial stiff- ness causes increase in central systolic and pulse pressures, which is a risk factor of cardiovascular events [1,137]. A number of physiological parame- ters have been proposed as indices of arterial stiffness, e.g. PWV, arterial distensibility and total arterial compliance. Of all these parameters, PWV is probably the most accepted index of arterial stiffness [123]. According to the Moens–Korteweg equation, PWV is related to the Young’s modulus (E) of a thin-walled elastic tube by: PWV = √ E× h/2rρ, where ρ is the fluid density, h is the wall thickness and r is the tube radius. [119,123]. It is generally accepted that PWV increases as arterial stiffness increases. By definition, PWV is the distance travelled by a (pressure or flow) pulse divided by the time of travel. Previous researchers have devised sev- eral methods to measure PWV by measuring PTT and travel distance. In this thesis, these previous methods are referred to as the direct methods as compared to the model-based method developed in this chapter. The direct methods use echocardiography [22,140], magnetic resonance imaging (MRI) [58,141] or applanation tonometry [112,130]. Using echocar- diography, PTT is measured as the delay between the flow pulses at the ascending and descending aorta while distance is measured on the long-axis image of the aortic arch. Using MRI, PTT is measured as the delay between the flow pulses at different locations along the entire aorta while distance is measured on the long-axis image of the aorta. Using applanation tonometry, PTT is measured as the delay between the pressure pulses at the carotid and femoral arteries while distance is measured on the surface of the skin. In the measurement of the carotid-femoral PWV using applanation tonom- etry, PTT is in the range of 30 to 100 ms while the distance is in the range of 40 to 100 cm [112,130]. In these ranges, a difference of 10 ms in PTT can re- sult in a difference of 100 cm/s in PWV. Therefore the precision requirement for PTT is rather high. A model-based method to measure PWV without PTT and distance measurements is developed in this chapter. A single uniform transmission line terminated by a complex load is used to model the arterial system (Figure 2.1). The applicability of this model has been described elsewhere [14,15]. The analogy between this model and the arterial system is explained in Appendix A. The line is assumed to be lossless since the resistance along the pulse transmission path is very small [15,25]. For a lossless transmission line, the phase velocity can be calculated from the line capacitance and the 18 2.1. Introduction C L R L R T L, C, l Figure 2.1: A transmission line terminated by a complex load: L is the line inductance; C is the line capacitance; l is the line length; CL is the load capacitance; RL is the load resistance; RT is the terminating resistance. line inductance. The phase velocity of the line corresponds to the PWV of the arterial system. To identify the model parameters, ascending aortic pressure (AAP) and ascending aortic flow (AAF) are used as the input and output to the model respectively. PWV is then calculated from the identified line capacitance and line inductance. The carotid artery pressure waveform is acquired non- invasively using applanation tonometry and used as a surrogate to the AAP waveform [81]. The AAF waveform is acquired using echocardiography. The purpose of this model-based method is to provide a tool to differ- entiate subjects with abnormal arterial stiffness from the healthy subjects. If this method can indicate abnormal arterial stiffness, it will be useful in predicting cardiovascular events and quantifying cardiovascular risks. In this thesis, a clinical study was performed to investigate the usefulness of this new method. In this study, the model-based method was applied to 8 healthy children and 8 children with Marfan syndrome (MFS). MFS is a connective tissue disorder. It is generally agreed that MFS causes abnormal elastic fiber content in the arterial wall, especially at the ascending aorta. Its primary cardiovascular manifestation is progressive dilatation of the aortic root eventually leading to aortic dissection or rupture [38]. It has been found in several studies that aortic pulse wave velocity is higher in subjects with MFS [22,58,59,140]. MFS is chosen in this study since its major manifestation is abnormal arterial stiffness. The objective was to investigate if the model-based method can differentiate subjects with MFS from healthy subjects. The model-based method was also compared to the direct method using echocardiography and the correlation between these two methods were found. 19 2.2. Methods 2.2 Methods 2.2.1 Outline of the model parameter identification algorithm The model has six parameters: line inductance (L), line capacitance (C), load capacitance (CL), load resistance (RL), terminating resistance (RT ) and line length (l). RL is dependent on L, C and RT , therefore there are only five free parameters. L and RT are estimated using some known relations. The remaining three free parameters, C, CL and l are found by minimizing a performance index (J). The carotid artery pressure waveform is used as the input to the model while the AAF waveform is the output. J consists of the difference between the model-predicted and the experimentally-measured AAF waveforms. An iterative optimization algorithm is used to minimize J . The transmission line model is numerically implemented using the finite element method as described in Appendix B. Since J may not be convex, a set of initial values for C, CL and l are chosen based on some physiological information of the subject. It is assumed that this set of initial values would lead to a local minimum that corresponds to the true values of C, CL and l. A version of the Newton conjugate-gradient method is used to minimize J [120]. It is a line search method, in which a direction is chosen from the current iterate and along this direction the next iterate is searched for such that J has a lower value. In this method, the direction is either the pure Newton direction or an approximation to it, found using the conjugate gradient method. Obtaining the Newton direction is a non-trivial task, and it can only be used as a search direction if the Hessian of J is positive definite. In the conjugate gradient method, the Newton direction is solved for iteratively. If the Hessian is found to be non-positive definite during the computation of the next iterate, the current iterate is chosen as an approximation to the Newton direction instead. The conjugate gradient method is very fast and therefore is chosen to compute the search direction. The simplest line search method is the steepest descent method, in which the search direction is the negative gradient of J at the current iterate. Although it is conceptually simple and easy to implement, the steepest descent method is often excruciatingly slow. In contrast, the Newton method typically has superlinear rate of convergence and is therefore chosen as the optimization algorithm in this study. 20 2.2. Methods 2.2.2 Estimations of the line inductance and the terminating resistance The transmission line equations are analogous to the one dimensional Navier- Stokes equations for incompressible, pulsatile flow in an elastic tube. The two set of equations are presented and the analogy is explained in Appendix A. The line inductance corresponds to the density of the fluid divided by the cross-sectional area of the tube. Accordingly, L is estimated as: L = ρb Aao (2.1) where ρb is the blood density andAao is the lumen area of the ascending aorta as measured by echocardiography. The standard value of blood density, 1.056 g/cm3 (or 0.00079 mmHgs2/cm2), is used. The terminating resistance corresponds to the peripheral resistance. Ac- cordingly, RT is estimated as mean (time-averaged) pressure (P̄ ) divided by mean (time-averaged) flow (Q̄): RT = P̄ Q̄ (2.2) 2.2.3 Initial values of the line capacitance, the load capacitance and the line length The line length corresponds to the distance from the root of the ascending aorta to the effective reflection site [27], which is also referred to as the effective arterial length. It has been suggested that the effective arterial length is dependent on age [114,119] and the relative distribution of the compliance in the arterial system [114]. Therefore, it is expected that l varies among different age groups and groups with different pathological conditions. l should also be dependent on body size since as height increases, the distances to each actual reflection sites also increase. In a previous invasive study [116], it was shown that the average effective arterial length was 48 cm for a group of normal men with the average body surface area (BSA) of 1.98 m2. In this method, a crude initial value is obtained by assuming that l varies linearly with BSA: linitial = BSA 1.98 × 48 (2.3) It should be noted that this equation was conceived to provide an initial value of l. In fact, l is dependent on a lot of factors and may not have a linear relationship with BSA. 21 2.2. Methods In the transmission line model, the total capacitance (CT ) is CT = C × l + CL (2.4) In this method, CT is initially estimated by the total arterial compliance (TAC) calculated using the area method [98], which is given by TAC = A12 RT (P1 − P2) (2.5) where P1 is the pressure at the very beginning of the diastolic decay, P2 is the end-diastolic pressure and A12 is the area under the pressure curve between P1 and P2 (Figure 2.2). In a previous modelling study using data collected from open-chest dogs [28], it was found that the total line capacitance (line capacitance times line length) was 40-50% of the total capacitance. Accord- ingly, the initial values of C and CL are calculated as follows: Cinitial = TAC 2linitial (2.6) C L,initial = TAC 2 (2.7) 2.2.4 Estimation of the load resistance At high frequencies, as viewed from the ventricles, the arterial system ap- pears to be minimally reflective [160]. It can be assumed that the termi- nating impedance matches the characteristic impedance of the line at high frequencies [28]. For a lossless line, the characteristic impedance is √ L/C and therefore RL is estimated as: RL = 1√ C L − 1 RT (2.8) 2.2.5 Performance index The performance index J is defined as follow: J = ∑n i (Qao,m(i)−Qao,e(i))2 n + µΛ (2.9) with Λ = (max (F −mn, 0))2 (2.10) 22 2.2. Methods 0 0.2 0.4 0.6 0.8 1 0 20 40 60 80 100 Time (s) P re ss u re ( m m H g ) P 1 P 2 A 12 Figure 2.2: The calculation of the total arterial compliance: P1 is the pres- sure at the beginning of the diastolic decay; P2 is the end-diastolic pressure; A12 is the area under the pressure curve between P1 and P2. 23 2.3. Experiments mn = 10 (2.11) where Qao,m is the model-predicted ascending aortic positive flow, Qao,e is the experimentally-measured ascending aortic positive flow, n is the number of data points, F is the magnitude of the model-predicted negative flow, mn is the allowed magnitude of the negative flow , and µ is the penalty coeffi- cient. In this study, the sampling rate for the pressure and flow waveforms was 12.0 kHz and n was about 104. The meanings of Qao,m and F are also displayed in Figure 2.3. In this method, when obtaining the experimentally-measured AAF from the Doppler spectrum, only the positive flow during the ventricular ejection is considered. The first term of J corresponds to mean-squared error be- tween the experimentally-measured and the model-predicted ascending aor- tic positive flow (Figure 2.3). A small amount of negative flow often occurs after the ventricular ejection. Therefore, J is defined to allow only a small amount of negative flow. The second term of J penalizes large magnitude of negative flow. The mean-squared error between the experimentally-measured and the model-predicted ascending aortic positive flows is expected to range from 102 to 1002 cm6/s. Accordingly, µ is chosen to be 1. If the magnitude of the negative flow is larger than 20, the penalty term would be of the same order of magnitude as the mean-squared error. Setting µ at 1 effectively restricts the search of the solution in a neighborhood with F smaller than 20 cm3/s. 2.2.6 Calculation of the model-based pulse wave velocity For a lossless transmission line, PWV is given by PWVmodel = 1√ LC (2.12) 2.3 Experiments 2.3.1 Data acquisition and processing This study was approved by the University of British Columbia Children’s and Women’s Research Ethics Board. 8 healthy children and 8 children with MFS participated in this study. 24 2.3. Experiments 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 −100 −50 0 50 100 150 200 250 300 350 Time (s) F lo w ( cm 3 /s ) Q ao,m |F| Figure 2.3: The model-predicted ascending aortic flow and its components: Qao,m is the positive ascending aortic flow; F is magnitude of the negative flow immediately after the ventricular ejection. Model-based pulse wave velocity Ascending aortic blood velocity was obtained using Doppler ultrasound with the transducer positioned at the suprasternal notch. The carotid artery pres- sure waveform was obtained simultaneously using applanation tonometry. Brachial artery systolic and diastolic blood pressures were measured using a cuff sphygmomanometer. The aortic annulus was measured using two- dimensional echocardiography. All the applanation-tonometric and echocar- diographic data were acquired using a GE Vivid 7 Pro ultrasound machine (GE Vingmed, Horten, Norway) equipped with an applanation tonometer. The AAF waveform was obtained by tracing the envelope of the Doppler spectrum of the ascending aortic blood velocity and multiplying it by the lu- men area (Figure 2.4, upper left panel). The lumen area was assumed to be circular and was calculated from the aortic annulus. Only the positive flow during the ventricular ejection was considered, and zero flow was assumed after the ventricular ejection (Figure 2.4, upper right panel). The carotid artery pressure waveform was used as a surrogate to the AAP waveform (Fig- ure 2.4, lower left panel). It has been previously shown that aortic pressure waveform measured invasively and carotid artery pressure waveform mea- 25 2.3. Experiments 0 0.5 1.0 1.5 2.0 0 0.5 1.0 1.5 2.0 −200 −100 0 100 200 300 Time (s) P re ss u re F lo w ( cm 3 /s ) 0 0.2 0.4 60 70 80 90 100 110 P re ss u re ( m m H g ) −100 0 100 200 300 400 F lo w ( cm 3 /s ) Time (s) Figure 2.4: Upper left panel: the Doppler spectrum and its envelope. Upper right panel: the averaged flow waveform; zero flow was assumed after the ventricular ejection. Lower left panel: the uncalibrated pressure waveform obtained using applanation tonometry. Lower right panel: the averaged pressure waveform. sured by applanation tonometry simultaneously have similar morphology [81]. The carotid artery pressure waveform was then calibrated by assuming that mean and diastolic pressures were constant throughout the arterial tree [82]. Mean pressure (MP) was calculated from the brachial artery systolic (SP) and diastolic pressures (DP) using the traditional formula [21]: MP = (SP+2DP)/3. Mean pressure was assigned to the area integral over the duration of the carotid artery pressure waveform, and diastolic pressure was assigned to the minimum of the carotid artery pressure waveform (Figure 2.4, lower right panel). Two sets of 4 to 6 beats were selected based on quality from each subject and averaged. The first set was used for parameter identification and the second set for model validation. 26 2.4. Results Direct pulse wave velocity using echocardiography The data acquisition procedure for the direct measurement of the aortic PWV using echocardiography has been described in details by Sandor et al. [140]. The Doppler spectrum of the ascending aortic blood velocity was acquired with the sample volume placed at the valve leaflets. Then, the Doppler spectrum of the descending aortic blood velocity was acquired with the sample volume placed as distally as possible. PTT was calculated as the time from the QRS to the onset of the descending aortic flow minus the time from the QRS to the onset of the AAF. The distance between the two recording sites was measured on the B-mode image of the aortic arch. 2.4 Results Data are expressed as mean±SD. All the clinical characteristics and the resulting estimates of the two groups were compared by means of an unpaired Student’s t-test. A value of p < 0.01 was considered statistically significant. The clinical characteristics of the subjects are listed in Table 2.1. TAC, the fixed model parameters including L and RT , and the initial model parameters including Cinitial, CL,initial and linitial, are listed in Table 2.2. Cinitial, CL,initial and linitial were calculated from TAC and BSA as describe in the previous section. The final model parameters, including C, CL and l are listed in Table 2.3. C, CL and l were obtained through the minimization of J . The ratio of C to CL is also listed in Table 2.3. The penalty coefficient (µ) is chosen to be 1 as described in the previous section. However, to investigate the effect of µ, the optimization of J was also performed with µ set at 0.1, 10 and 100. The resulting PWVµ=0.1, PWVµ=1, PWVµ=10, PWVµ=100 are listed in Table 2.4. PWVecho, which is the direct PWV using echocardiography, and PWVinitial, which is calcu- lated using Cinitial and L are also listed in Table 2.4. Model validation are evaluated by the root-mean-squared error (RMSE) between the model-predicted and experimentally-measured ascending aortic positive flow waveforms. Two RMSE values were calculated with one us- ing Cinitial, CL,initial and linitial, and one using C, CL and l. The final estimates were found by setting µ at 1. These values are listed in Table 2.6. A sample of the model-predicted AAF waveform and the corresponding experimentally-measured AAF waveform is shown in Figure 2.8. The box and whisker plots of PWVµ=1 and PWVecho are shown in Figure 2.5. The correlation between PWVµ=1 and PWVecho is shown in 27 2.4. Results Control (8) MFS (8) p-value Age (yrs) 11± 2 15± 3 NS BSA (m2) 1.1± 0.2 1.7 ± 0.2 p < 0.0006 SPBA (mmHg) 107 ± 8 111 ± 7 NS DPBA (mmHg) 63± 5 70± 7 NS AO (cm) 1.62 ± 0.21 2.16 ± 0.25 p < 0.0003 Peak AAF (cm3/s) 265 ± 74 328 ± 117 NS Table 2.1: Clinical characteristics of the subjects: BSA is the body surface area; SPBA is the systolic pressure at the brachial artery; DPBA is the diastolic pressure at the brachial artery; AO is the aortic annulus; Peak AAF is the peak ascending aortic flow. Control (8) MFS (8) p-value TAC (cm3/mmHg) 1.2± 0.4 2.3± 1.0 NS L (mmHgs2/cm3/cm) 0.00041 ± 0.00012 0.00022 ± 0.00007 p < 0.003 RT (mmHgs/cm 3) 1.3± 0.3 1.3± 0.6 NS Cinitial (cm 3/mmHg/cm) 0.022 ± 0.007 0.028 ± 0.009 NS C L,initial (cm 3/mmHg) 0.61 ± 0.19 1.14 ± 0.50 NS linitial (cm) 28± 6 40± 5 p < 0.0006 Table 2.2: Fixed model parameters and initial model parameters: TAC is the total arterial compliance; L is the line inductance; RT is the terminating resistance; Cinitial is the initial value of the line capacitance; CL,initial is the initial value of the load capacitance; linitial is the initial value of the line length. Figure 2.6. The Bland-Altman plot of PWVµ=1 and PWVecho is shown in Figure 2.7. The model-based PWV with different µ were found to be significantly higher in the MFS group. Both PWVinitial and PWVecho are not signifi- cantly different between the MFS and the control groups. The correlations between PWVecho and different model-based PWV are shown in Table 2.5 28 2.4. Results Control (8) MFS (8) p-value C (cm3/mmHg/cm) 0.015 ± 0.007 0.013 ± 0.006 NS CL (cm 3/mmHg) 0.63 ± 0.19 1.35 ± 0.69 NS l (cm) 25± 7 34± 5 NS C/CL 0.025 ± 0.009 0.011 ± 0.003 p < 0.003 Table 2.3: Final model parameters: C is the final estimate of the line capac- itance; CL is the final estimate of the load capacitance; l is the final estimate of the line length. Control (8) MFS (8) p-value PWVµ=1 (cm/s) 438 ± 97 626 ± 115 p < 0.004 PWVµ=0.1 (cm/s) 383 ± 112 578 ± 143 p < 0.01 PWVµ=10 (cm/s) 410 ± 92 600 ± 119 p < 0.004 PWVµ=100 (cm/s) 398 ± 107 605 ± 132 p < 0.005 PWVecho (cm/s) 362 ± 133 444 ± 115 NS PWVinitial (cm/s) 349 ± 69 422 ± 55 NS Table 2.4: Pulse wave velocity: PWVmu=1 is the model-based pulse wave ve- locity with the penalty coefficient set at 1; PWVmu=0.1 is the model-based pulse wave velocity with the penalty coefficient set at 0.1; PWVmu=10 is the model-based pulse wave velocity with the penalty coefficient set at 10; PWVmu=100 is the model-based pulse wave velocity with the penalty coeffi- cient set at 100; PWVecho is the direct measurement of pulse wave velocity using echocardiography; PWVinitial is the pulse wave velocity estimated using the initial value of the line capacitance. 29 2.4. Results Correlation PWVecho vs PWVµ=1 0.51 PWVecho vs PWVµ=0.1 0.26 PWVecho vs PWVµ=10 0.60 PWVecho vs PWVµ=100 0.52 PWVecho vs PWVinitial 0.50 Table 2.5: Correlations between the echocardiographic pulse wave velocity and different model-based pulse wave velocities: PWVmu=1 is the model- based pulse wave velocity with the penalty coefficient set at 1; PWVmu=0.1 is the model-based pulse wave velocity with the penalty coefficient set at 0.1; PWVmu=10 is the model-based pulse wave velocity with the penalty coefficient set at 10; PWVmu=100 is the model-based pulse wave velocity with the penalty coefficient set at 100; PWVecho is the direct measurement of pulse wave velocity using echocardiography; PWVinitial is the pulse wave velocity estimated using the initial value of the line capacitance. Control (8) MFS (8) RMSEinitial (cm 3/s) 45.2 ± 18.4 92.4 ± 40.7 RMSE (cm3/s) 27.4 ± 6.2 31.0± 17.30 Table 2.6: Performance of the model: RMSE is the root-mean-squared error between the model-predicted and experimentally-measured ascending aortic positive flow waveforms; RMSEinitial was calculated using the initial values of the line capacitance, load capacitance and line length; RMSE was cal- culated using the final estimates of the line capacitance, load capacitance and line length. The final estimates were found using optimization with the penalty coefficient set at 1. 30 2.4. Results 200 250 300 350 400 450 500 550 600 650 700 ( cm /s ) PWV μ=1 Control PWV echo Control PWV μ=1 MFS PWV echo MFS model-based PWV echocardiographic PWV Figure 2.5: Box and whisker plots: medians, upper and lower quantiles and ranges are displayed; PWVµ=1 is the model-based pulse wave velocity; PWVecho is the direct pulse wave velocity obtained using echocardiography. 31 2.4. Results 0 200 400 600 800 1000 0 100 200 300 400 500 600 700 800 900 1000 P W V μ = 1 (c m /s ) PWVecho (cm/s) r = 0.51 y = 0.57x+303 Figure 2.6: Correlation between the model-based pulse wave velocity and the direct pulse wave velocity: PWVµ=1 is the model-based pulse wave velocity; PWVecho is the direct pulse wave velocity obtained using echocardiography. 32 2.4. Results 200 300 400 500 600 700 −400 −300 −200 −100 0 100 200 300 400 500 600 mean + 1.96 SD mean - 1.96 SD mean P W V μ = 1 P W V ec h o ( cm /s ) (PWV μ=1 + PWV echo )/2 (cm/s) Figure 2.7: Bland-Altman plot of the model-based pulse wave velocity and the direct pulse wave velocity: PWVµ=1 is the model-based pulse wave velocity; PWVecho is the direct pulse wave velocity obtained using echocar- diography. 33 2.4. Results 0 0.2 0.4 0.6 −100 0 100 200 300 400 F lo w ( cm 3 /s ) Time (s) Model-predicted Experimentally-measured Figure 2.8: A sample of the model-predicted ascending aortic flow wave- form and the corresponding experimentally-measured ascending aortic flow waveform: the root-mean-squared error between the model-predicted and experimentally-measured ascending positive flow waveforms is 16.6 cm3/s. 34 2.5. Discussion 2.5 Discussion 2.5.1 Ability to differentiate the subjects with Marfan syndrome from the healthy subjects PWVµ=1, PWVµ=0.1, PWVµ=10 and PWVµ=100 were higher in the MFS group. It suggests that the proposed model-based pulse wave velocity is able to differentiate the children with MFS from the normal children regardless of the value of the penalty coefficient. PWVinitial was not significantly different between the two groups. It suggests that the optimization step is necessary to obtain a PWV estimate that can differentiate the children with MFS from the normal children. PWVecho was also not significantly different between the two groups. This is unexpected since pulse wave velocity obtained using echocardiog- raphy was found to be higher in subjects with MFS in two other studies [22,140]. A possible reason for this is the small sample size, and therefore the power of the test was small. In the box and whisker plots (Figure 2.5), it is clear that PWVµ=1 of the control group and PWVµ=1 of the MFS group occupied different ranges. On the other hand, for PWVecho, the range for the control group and the range for the MFS group largely overlapped. It is also clear that, for the MFS group, PWVµ=1 is higher than PWVecho. A possible reason for this is the larger aortic diameter in the MFS group (Table 2.1). By equation (2.1), larger aortic diameter gives smaller L. And by equation (2.12), smaller L gives larger PWV. However, it should not be concluded that the ability of PWVµ=1 to differentiate the two groups was solely due to L. It should be noted that PWVinitial has the same dependency on the aortic diameter but it lacks the ability to differentiate the two groups. To determine if PWVµ=1 can detect abnormal arterial stiffness indepen- dent of the size of aorta, a study on two groups with the aortic diameter matched is necessary. 2.5.2 Correlations between the pulse wave velocities obtained by different methods The correlation coefficient between PWVµ=1 and PWVecho was 0.51 (Figure 2.6). To put this into perspective, a previous study that compared PWV obtained by three devices: Complior, SphygmoCor, and Arteriograph is considered [130]. Both Complior and SphygmoCor measure carotid-femoral PWV, but the algorithms to determine PTT and distance are different. The 35 2.5. Discussion main principle behind Arteriograph is to occlude the brachial artery by ap- plying a cuff pressure exceeding systolic pressure by 35–40 mmHg. PTT is then measured as the time interval between the incident and reflected sys- tolic wave. The correlation coefficients between Complior and SphygmoCor, Complior and Arteriograph, SphygmoCor and Arteriograph were 0.57, 0.36, 0.29 respectively. The correlation between the proposed model-based PWV and the echocardiographic PWV is similar to that between Complior and SphygmoCor. 2.5.3 Interchangeability between the pulse wave velocities obtained by different methods The limits of agreement between PWVµ=1 and PWVecho were -133 cm/s and 391 cm/s. Therefore, PWVµ=1 and PWVecho are not interchangeable. Similar conclusion for Complior and SphygmoCor was found by Millasseau et al. [112]. Both Complior and SphygmoCor are validated techniques to estimate carotid-femoral PWV, and yet they are not interchangeable. It suggests that lack of interchangeability does not affect the validity of either technique. In fact, the PWV estimates are strongly dependent on the underlying algorithms [112,130]. It is recommended that the same PWV estimation method is used when differentiating two different populations. 2.5.4 Performance of the model in predicting ascending aortic flow The experimentally-measured peak AAF was around 300 cm3/s (Table 2.1). Relative to the peak AAF, the RMSE was reasonable (Table 2.6). As shown in Figure 2.8, the model-predicted AAF waveform is reasonably close to the experimentally-measured AAF waveform. In addition, using the initial values, the model-predicted ascending aor- tic positive flow waveform had poor agreement with the experimentally- measured ascending aortic positive flow waveform (Table 2.6). It suggests that the optimization step is necessary. 2.5.5 Large artery compliance to small artery compliance ratio It is interesting to note that the ratio of C to CL was lower in the MFS group. In the transmission line model, the line can be considered as representing the aorta and the elastic arteries with the load representing the muscular 36 2.6. Conclusions arteries and the arterioles. It has been previously shown that increase in arterial stiffness is largely confined to the aorta in the MFS patients [76]. Therefore, it seems reasonable to expect that the ratio of C to CL would be lower with the MFS patients. Indeed, the ratio of C to CL found in this study matched this expectation. It suggests that this model and the corresponding identification method can also model the difference in the distribution of compliance. 2.6 Conclusions In this chapter, a model-based method to estimate PWV is presented. This method was compared to the direct measurement of PWV using echocar- diography. A few conclusions were made. First, the model-based PWV had stronger ability than the direct PWV obtained using echocardiography in differentiating the children with MFS from the normal children. Second, the model-based PWV and the direct PWV obtained using echocardiogra- phy showed a moderate correlation. Third, the model-based PWV and the direct PWV obtained using echocardiography are not interchangeable. The number of subjects in this study was small. A follow-up study including a larger number of subjects is necessary to confirm the clinical usefulness of this model-based method. 37 Chapter 3 Monitoring Cardiac Output using Grey-Box Identification on a State-Space Representation of the Transmission Line Model1 Cardiac output (CO) is a strong indicator of the adequacy of the cardio- vascular circulation. CO normalized to the body surface area is associated with the survival rate of postoperative patients [18,40]. A minimally-invasive technique that allows continuous monitoring of CO is valuable. It has been suggested that CO can be estimated from the periph- eral pressure waveform using the pulse contour analysis. In the traditional pulse contour techniques, CO is related to the area of the systolic portion of the pressure waveform. In this chapter, a novel pulse contour technique that uses the periph- eral pressure waveform and a transmission line model to estimate CO is developed. The transmission line model represents the pulse transmission path from the ascending aorta to the radial artery. Calibration is performed through the acquisition of the ascending aortic flow waveform using echocar- diography. A patient-specific transmission line model is identified from the ascending aortic flow and the radial artery pressure. With the identified model and the radial artery pressure waveform, the ascending aortic flow is estimated. CO is calculated from the estimated ascending aortic flow. 1 A version of this chapter has been submitted. 38 3.1. Introduction 3.1 Introduction CO measurement during critical care is traditionally performed using ther- modilution, the direct Fick method or dye-dilution. All these methods re- quire the placement of a pulmonary artery catheter, which is associated with increased risk and increased utilization of resources [34]. A less inva- sive alternative is much desired. CO can also be measured non-invasively using transthoracic echocardiography [144], [145], [169]. However, a trained technician is required to perform the measurement. Therefore, continuous CO monitoring in the intensive care unit (ICU) is not commonly performed using transthoracic echocardiography. In ICU, patients often have an indwelling radial artery catheter inserted to monitor the vital signs. Radial artery pressure (RAP) waveform is there- fore readily available. For many years, researchers have tried to derive CO from the peripheral artery pressure (PAP) waveform. Various pressure wave- form analysis techniques, which are also referred to as the pulse contour anal- ysis techniques, have been put forward [54], [74], [97], [115], [173]. Many of these techniques involve one or several measurements of CO before the monitoring commences, and a relationship between the CO values and the PAP waveform is thus derived. Such calibration is often performed using transpulmonary lithium-dilution. In this chapter, a new CO estimation technique is developed. According to this new technique, the pulse transmission path from the ascending aorta to the radial artery is represented by a transmission line model. During calibration, the ascending aortic flow (AAF) waveform is acquired using transthoracic echocardiography. The transmission line model parameters are identified using the AAF and the RAP waveforms. The new CO estimation technique is outlined as follow: 1. Calibration (a) An empirical state-space (SS) model is identified from the AAF and the RAP waveforms using subspace identification. (b) The transmission line model parameters are then derived from the empirical SS model. 2. Estimation (a) The transfer function from the AAF to the RAP is calculated from the identified model parameters. 39 3.1. Introduction C L R L R T L, C, l Figure 3.1: A lossless transmission line terminated by a complex load: L is the line inductance; C is the line capacitance; l is the line length; CL is the load capacitance; RL is the load resistance; RT is the terminating resistance. (b) The AAF is estimated from the RAP through the inverse of the identified transfer function (c) CO is calculated from the area-time integral of the model-estimated AAF waveform. This new technique has two advantages compared to previous techniques. First, the calibration is performed using transthoracic echocardiography and is completely non-invasive. Second, assuming that the flow waveform con- tains more information than the CO value, then the identified model should be able to describe the pressure-flow dynamics more precisely. The transmission line model consists of a single uniform lossless line ter- minated by a complex load (Fig. 3.1). The applicability of this model has been described elsewhere [14], [15]. The analogy between this model and the pulse transmission path is explained in Appendix A. In this chapter, the electric terms are used instead of the fluid dynamic terms. A SS representa- tion of this model based on the finite-element method [96] is developed. In this representation, the unknown entries of the system matrices are linear functions of (the reciprocals of) the transmission line model parameters. An empirical SS model can be identified from the input and output data using subspace identification [39], [99], but it is often more useful to identify a SS model with physical meaning. The reason for this is one may be inter- ested in identifying the physical parameters, or one may want to incorporate additional knowledge about the system into the identification process. Re- cently, Kind et al. [85] used the three- and four-element Windkessel model to describe the arterial system, and proposed a method to derive the Wind- kessel parameters from the empirical SS model. This method was specifically designed for the three- and four-element Windkessel model. Xie and Ljung [175] proposed a more general framework to identify a physical SS model 40 3.1. Introduction from an empirical SS model. In this framework, the transformation matrix between the empirical SS model and the physical SS model is found through iterative bilinear regression. It was suggested that this framework is suitable when the system matrices of the physical SS model have a small number of unknown entries or a special structure [102], [175]. In this chapter, it is proposed that this framework is also suitable when the unknown entries of the system matrices are linear functions of a small number of the physical model parameters. An identification method based on this framework is developed. Since this method involves black-box identi- fication (subspace identification) and a model structure with physical mean- ing (the transmission line model), it is a grey-box identification method. In the Methods section, the followings are presented: 1. the derivation of the SS representation of the transmission line model, 2. the derivation of the physiological ranges of the physical model pa- rameters, 3. the local identifiability of the proposed SS model, 4. the grey-box identification method, and 5. the estimated bias and variances of the proposed grey-box identifica- tion method, 6. the performance of the proposed grey-box identification method on simulated data, 7. the CO estimation technique. In the Results section, the followings are presented: • the limits of agreement between the proposed CO estimation technique and the more well-established CO estimation technique using transtho- racic echocardiography. The data were collected from 5 infant patients during their stay in ICU after cardiac surgeries. 41 3.2. Methods 3.2 Methods 3.2.1 Derivation of the state-space representation of the transmission line model Equations of a lossless transmission line terminated by a complex load The voltage and the current of a transmission line are described by the following equations: − ∂V (x, t) ∂x = L ∂I(x, t) ∂t (3.1) −∂I(x, t) ∂x = C ∂V (x, t) ∂t (3.2) where V is the voltage, I is the current, L is the line inductance, C is the line capacitance, x is the spatial index and t is the time index. The transmission line equations are analogous to the one-dimensional Navier-Stokes equations with the voltage corresponds to the pressure and the current corresponds to the flow (Appendix A). The complex load at the end of the line is described by the following equations: CL dΦc dt = I(l, t)− V (l, t) RT (3.3) I(l, t) = ( 1 RT + 1 RL ) V (l, t)− 1 RL Φc (3.4) where CL is the load capacitance, RL is the load resistance, RT is the ter- minating resistance, l is the length of the line and Φc is the voltage across CL. Finite element representation of the transmission line model In the finite element method, voltage and current are evaluated at a finite number of nodes. Each voltage node is placed between two current nodes. Voltage is evaluated at half a time step before and after current. The trans- mission line equations (3.1) and (3.2) can be approximated using the finite element method as described by Lee et al. [96]: V n+1j = K3V n j −K4(I n+ 1 2 j+1 − I n+ 1 2 j ) (3.5) I n+ 1 2 j+1 = K1I n− 1 2 j+1 −K2(V nj+1 − V nj ) (3.6) 42 3.2. Methods Spatial nodes I 1 I 2 V 1 V 2 ... ... V j V j+1 I j I j+1 ... V p-1 I p-1 I p ... Time steps I 1/2 I 1+1/2 V 1 V 2 ... ... V n V n+1 I n-1/2 I n+1/2 ... ... Time ∆x ∆x ∆t Figure 3.2: The indexing for the current and the voltage in the finite element representation: the superscript indicates the time step and the subscript indicates the spatial node; n is used to indicates the time step; j is used to indicates the spatial node; p − 1 indicates the last voltage node and p indicates the last current node; ∆x is the distance between two current nodes or two voltage nodes; ∆t is the duration of each time step. where j indicates the node and n indicate the time step (Fig. 3.2). K1 to K4 are given by K1 = 1 (3.7) K2 = ∆t L∆x (3.8) K3 = 1 (3.9) K4 = ∆t C∆x (3.10) For numerically stability, this approximation scheme requires [96]: µ ≤ ∆x ∆t (3.11) where µ is the wave propagation speed. For a lossless transmission line, µ is given by: µ = 1 LC (3.12) 43 3.2. Methods Equations (3.3) and (3.4) are approximated as follows: Φc n+1 = ∆t CL I n+ 1 2 p − ∆t 2CLRT ( V n+1p−1 + V n p−1 ) + Φc n (3.13) I n+ 1 2 p = ( 1 RT + 1 RL ) V np−1 − 1 RL Φc n (3.14) where p− 1 indicates the last voltage node and p indicates the last current node (Fig. 3.2). State-space representation of the finite element transmission line model Since the time index for I is arbitrary, In− 1 2 and In+ 1 2 can be re-indexed as In and In+1 respectively. Equations (3.5), (3.6), (3.13) and (3.14) are rewritten as follow: V n+1j = K3V n i −K4(In+1j+1 − In+1j ) (3.15) In+1j+1 = K1I n j+1 −K2(V nj+1 − V nj ) (3.16) Φc n+1 = ∆t CL In+1p − ∆t 2CLRT ( V n+1p−1 + V n p−1 ) + Φc n (3.17) In+1p = ( 1 RT + 1 RL ) V np−1 − 1 RL Φc n (3.18) It should be noted that In represents current that is evaluated at half a time step later than V n. By approximating, on the right-hand side, In+1 with In, V n+1p−1 with V n p−1, and V np−1 with V n−1 p−1 , equations (3.15) and (3.17) can be written as V n+1j = K3V n j −K4(Inj+1 − Inj ) (3.19) Φc n+1 = ∆t CL Inp − ∆t 2CLRT ( V np−1 + V n−1 p−1 ) +Φc n (3.20) Intuitively, this approximation requires ∆t to be small. Similar requirement on ∆t is also imposed by the inequality in (3.11). Equations (3.16), (3.18)-(3.20), with the current at the beginning of the line as the input (u) and the voltage at the end of the line as the output 44 3.2. Methods (y), can then be written in the SS form: x(n+ 1) = Ax(n) +Bu(n) (3.21) y(n) = Cx(n) (3.22) with u(n) = I1(n) = I n 1 (3.23) y(n) = Vp−1(n) = V n p−1 (3.24) x(n) =   V1(n) I2(n) V2(n) I3(n) ... Vp−1(n) Ip(n) Φc(n) Vp−1(n − 1)   =   V n1 In2 V n2 In3 ... V np−1 Inp Φc n V n−1p−1   (3.25) A =   K3 −K4 0 · · · 0 K2 K1 −K2 0 · · · 0 0 K4 K3 −K4 0 · · · 0 ... . . . . . . . . . . . . . . . · · · ... 0 · · · 0 K2 K1 −K2 0 0 0 0 · · · 0 K4 K3 −K4 0 0 0 · · · 0 (α+A1) 0 −A1 0 0 · · · 0 −αA22 A2 1 −αA22 0 · · · 0 1 0 0 0   (3.26) B =   K4 0 ... 0   (3.27) C = [ 0 · · · 0 1 0 0 0 ] (3.28) 45 3.2. Methods where α = 1 RT (3.29) A1 = 1 RL (3.30) A2 = ∆t CL (3.31) The system order is 2p. The terminating resistance RT corresponds to the peripheral resistance. Peripheral resistance is calculated as the mean aortic pressure divided by the cardiac output [119]. It is often assumed that mean pressure stays the same along the arterial tree. Cardiac output is the same as mean aortic flow. Therefore, peripheral resistance can also be approximated as the mean pressure at any location along the arterial tree divided by the mean aortic flow. In analogy, RT can be approximated as: RT = V̄p−1 Ī1 (3.32) The remaining unknowns in the system matrices are: θ = [ K2 K4 A1 A2 ] (3.33) and the known entries are: ψ = [ K1 K3 1 ] (3.34) Since ∆x is l/(p − 1), the total line inductance (Ltot) and the total line capacitance (Ctot) can be retrieved from K2 and K4: Ltot = L× l = ∆t(p− 1) K2 (3.35) Ctot = C × l = ∆t(p− 1) K4 (3.36) RL and CL can be retrieved from A1 and A2. 3.2.2 Physiological ranges of the physical model parameters The physiological ranges of Ctot, CL, the characteristic impedance of the line (Zch) and the ratio Zch/RT are given in Table 3.1. These ranges were derived from the values found in literature. 46 3.2. Methods As viewed from the ventricles, the arterial system appears to be mini- mally reflective at high frequencies [160]. Therefore, it can be assumed that Zch matches the load impedance at high frequencies. The following relation thus holds: 1 RL = 1 Zch − 1 RT (3.37) Traditionally, Zch is estimated from the input impedance, as calculated from the aortic pressure and flow, of the arterial system [119]. It has been shown that Zch for a group of children aged 0.8 to 4 years was in the range of 0.08–0.6 mmHgs/cm3, and Zch decreases with age [26]. It is then deduced that Zch for infants is in the range of 0.08-2 mmHgs/cm 3. The total arterial compliance (TAC) defined by the Windkessel model represents the collective compliance of the arterial system [98] . In the trans- mission line model, the sum of Ctot and CL also represents the collective compliance of the arterial system. Therefore, it can be assumed that TAC approximates the sum of Ctot and CL. It has been shown TAC for a group of children aged 3 to 12 years was 0.71 cm3/mmHg [139]. It has also been shown that TAC increases with age among children because of the increase in the arterial size [143]. It is then deduced that TAC for infants is in the range of 0.01–0.7 cm3/mmHg. In a previous modelling study on data col- lected from open-chest dogs [28], it was found that the ratio of Ctot to CL was about 0.7 to 1. Assuming the same ratio is applicable to human infants, Ctot and CL for infants are in the range of 0.01–0.5 cm 3/mmHg. It has also been shown that the ratio of Zch to RT for a group of children aged 0.8 to 4 years was in the range of 0.05–0.23, and this ratio decreases with age [26]. It is then deduced that the Zch to RT ratio for infants is in the range 0.05–0.5. The ranges of Ltot, RT and RL can be calculated from the ranges of Ctot, Zch and Zch/RT . 3.2.3 Local identifiability of the proposed state-space model In the proposed SS model, the unknown physical model parameters are Ctot, Ltot, CL and RL. To investigate the local identifiability of this model with respect to these parameters, a method proposed by Dötsch and Van den Hof [42] was used. If the model described by equations (3.21) to (3.31) is stable, then the rank of CdPC T d (3.38) 47 3.2. Methods Ctot 0.01–0.5 cm 3/mmHg CL 0.01–0.5 cm 3/mmHg Zch 0.08–2 mmHgs/cm 3 Zch/RT 0.05–0.5 mmHgs/cm 3 Table 3.1: Physiological ranges of the physical model parameters: Ctot is the total line capacitance; CL is the load capacitance; Zch is the characteristic impedance; RT is the terminating resistance. equals the number of parameters if and only if the model is locally identi- fiable with respect to these parameters. P is the solution to the Lyapunov equation P = AdPA T d +BdB T d (3.39) where Ad, Bd and Cd are given by Ad =   A 0 0 0 0 ∂A ∂Ctot A 0 0 0 ∂A ∂Ltot 0 A 0 0 ∂A ∂CL 0 0 A 0 ∂A ∂RL 0 0 0 A   (3.40) Bd =   B ∂B ∂Ctot ∂B ∂Ltot ∂B ∂CL ∂B ∂RL   (3.41) Cd =   ∂C ∂Ctot C 0 0 0 ∂C ∂Ltot 0 C 0 0 ∂C ∂CL 0 0 C 0 ∂C ∂RL 0 0 0 C   (3.42) (3.43) For a 2-element transmission line model, 48 3.2. Methods ∂A ∂Ctot =   0 2∆t C2 tot 0 0 0 0 0 0 0 0 0 0 0 − 2∆t C2 tot 0 2∆t C2 tot 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0   (3.44) ∂A ∂Ltot =   0 0 0 0 0 0 − 2∆t L2 tot 0 2∆t L2 tot 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0   (3.45) ∂A ∂CL =   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ∆t 2αC2 L −∆t C2 L 0 ∆t 2αC2 L 0 0 0 0 0 0   (3.46) ∂A ∂RL =   0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 − 1 R2 L 0 1 R2 L 0 0 0 0 0 0 0 0 0 0 0 0 0   (3.47) ∂B ∂Ctot =   − 2∆t C2 tot 0 0 0 0 0   (3.48) 49 3.2. Methods ∂B ∂Ltot = ∂B ∂CL = ∂B ∂RL =   0 0 0 0 0 0   (3.49) ∂C ∂Ctot = ∂C ∂Ltot = ∂C ∂CL = ∂C ∂RL (3.50) = [ 0 0 0 0 0 0 ] (3.51) Ad, Bd and Cd for model with more elements can be derived easily. Local identifiability test for the 2-, 3-, 4- and 5-element models were performed. Values of Ctot, CL, Zch and Zch/RT were randomly chosen from the ranges given in Table 3.1. Ltot, RT and RL were calculated from the chosen values. The eigenvalues of A and the rank of CdPC T d were found. The model was stable if all the eigenvalues ofA were inside the unit circle. 5000 different sets of parameter values were tested. ∆t was chosen to be 10−4s, 10−5s and 10−6s. It was found that some of the models were unstable (Table 3.2). The number of unstable models increases as ∆t increases, and as number of the elements increases. This may be related to the requirement on ∆t and ∆x given by the inequality in (3.11). However, it is reasonable to assume that the pulse transmission path from the AAF to the RAP is a stable system and would be represented by a stable model. In all cases, if the model was stable, the rank of CdPC T d was 4. 3.2.4 Grey-box identification Framework of the grey-box identification The physical SS model (A(θ), B(θ), C(θ)) that describes a lossless transmis- sion line terminated by a complex load has been described in the previous section and is given by equations (3.21)-(3.31). An empirical SS model (Â, B̂, Ĉ) is identified from the input, I1(n), and the output, Vp−1(n), using subspace identification. The transformation matrix T relates the subspace- identified empirical SS model and the physical SS model as follows: T−1A(θ)T = Â (3.52) 50 3.2. Methods Number of stable model ∆t = 10−4, 2-element 5000 ∆t = 10−4, 3-element 4942 ∆t = 10−4, 4-element 4606 ∆t = 10−4, 5-element 3798 ∆t = 10−5, 2-element 5000 ∆t = 10−5, 3-element 5000 ∆t = 10−5, 4-element 4990 ∆t = 10−5, 5-element 4954 ∆t = 10−6, 2-element 5000 ∆t = 10−6, 3-element 5000 ∆t = 10−6, 4-element 5000 ∆t = 10−6, 5-element 5000 Table 3.2: The number of stable models for different time steps and numbers of elements: ∆t is the time step. T−1B(θ) = B̂ (3.53) C(θ)T = Ĉ (3.54) In the framework proposed by Xie and Ljung [175], T and θ are found by minimizing the error of equations (3.52)-(3.54) through iterative bilinear regression. The error criterion is defined as e = ∥∥∥TÂ−A(θ)T∥∥∥2 F + ∥∥∥B(θ)−TB̂∥∥∥2 F + ∥∥∥Ĉ−C(θ)T∥∥∥2 F (3.55) where ‖·‖ is the Frobenius norm. T and θ are updated alternatively by applying linear regression on equations (3.52)-(3.54). The iterative search begins with an initial estimate of θ. Subspace identification using N4SID The N4SID algorithm has been described in [165]. In brief, an extended observability matrix is determined from the projection of the future outputs onto the past and future inputs and the past outputs. The system matrices are then determined from a set of linear equations, which is written in terms of the extended observability matrix, inputs and outputs. The Matlab algorithm N4SID was used in this thesis. 51 3.2. Methods Derivation of θ from Â, B̂, Ĉ Equations (3.52)-(3.54) can be written in terms of the vectorized T as fol- lows:   Â T ⊗ I − I ⊗A(θ) B̂T ⊗ I I ⊗C(θ)   vec(T) =   0B(θ) vec(Ĉ)   (3.56) where ⊗ denotes the Kronecker product. Equations (3.52) and (3.53) can also be written in terms of the vectorized A(θ) and B(θ) respectively [TT ⊗ I]vec(A(θ)) = TÂ (3.57) B(θ) = TB̂ (3.58) Since the entries in A(θ) are linear functions of θ, two matrices HA and GA can be defined such that vec(A(θ)) = HAθ T +GAψ T (3.59) with vec(A(θ))(i) = HA(i,1)K2 +HA(i,2)K4 +HA(i,3)A1 +HA(i,4)A2 + GA(i,1)K1 +GA(i,2)K3 +GA(i,3) (3.60) where (i) indicates the element at the ith row, and (i, w) indicates the element at the ith row and wth column. Two matrices HB and GB are defined similarly: B(θ) = HBθ T +GBψ T (3.61) Equations (3.57) and (3.58) can then be written in terms of θ as follow[ TT ⊗ I 0 0 I ] [ HA HB ] θT = [ vec(TÂ) TB̂ ] − [ TT ⊗ I 0 0 I ] [ GA GB ] ψT (3.62) To minimize e given in equation (3.55), iterative bilinear regressions are performed as follows: 52 3.2. Methods 1. A set of initial values for θ are chosen. 2. T is solved for in equation (3.56) using linear regression. 3. With the updated T, θ is solved for in equation (3.62) using linear regression. 4. Step 2 and 3 are repeated alternatively with the updated θ and T until the gradient of e is smaller than some chosen bound. An additional nonlinear constraint Zch of a lossless line is given by Zch = √ L C (3.63) As viewed from the ventricles, the arterial system appears to be minimally reflective at high frequencies [160]. Therefore, it can be assumed that the characteristic impedance of the line matches the load impedance at high frequencies. The following relation thus holds: 1 RL = 1 Zch − 1 RT (3.64) A1 = √ K2 K4 − α (3.65) To implement this nonlinear constraint, HA, HB, GA andGB are modified. Let K−2 and K − 4 be the values of K2 and K4 found in the previous iteration. Then A1 is approximated as A1 = √ K−2 K−4 − α (3.66) A new θ is defined as θ = [ K2 K4 A2 ] (3.67) and a new ψ is defined as ψ = [ K1 K3 1 A1 ] (3.68) HA, HB, GA and GB are re-defined accordingly. To find T and θ, iterative bilinear regression is applied on equation(3.56) and equation (3.62) with the re-defined HA, HB, GA and GB. 53 3.2. Methods ˆ ˆ ˆ θ = [K 2 K 4 A 2 ] A(θ), B(θ), C(θ) are initialized using R T and the inital guess of θ H A, H B, G A, G B are defined using equations (3.44) and (3.46) Iterative bilinear regressions are applied on equations (3.41) and (3.47) to find θ A, B, C are identified from I 1 and V p-1 using N4SID subspace identification Figure 3.3: The grey-box identification algorithm: Â, B̂, Ĉ represent the empirical model; I1 is the current at the beginning of the line or analogously the ascending aortic flow; Vp−1 is the voltage at the end of the line or or analogously the radial artery pressure; K2, K4 and A2 are the unknowns of the physical model; A(θ), B(θ), C(θ) represent the physical model; RT is the terminating resistance. Summary of the grey-box identification algorithm The steps of the identification algorithm are shown in Figure 3.3. 3.2.5 Estimation of the bias and variance using a Monte-Carlo study A Monte-Carlo study was performed to estimate the bias and the variance of the proposed grey-box identification method. Several sets of parameter val- ues were used as shown in Table 3.3. ∆t was 10−4. RT was 10 mmHgs/cm 3. Number of elements was 3. The input was the realization of a Gaussian random variable with zero mean and unit variance. The length of the realization was 105 and 100 realizations were used. 54 3.2. Methods The estimates of the bias ( ˆBias) and the variance (V̂ar) were calculated as follows: ϕ = [ Ctot Ltot CL ] (3.69) ˆBias = 1 M ∣∣∣∣∣ M∑ i=1 (ϕ0 − ϕ̂i) ∣∣∣∣∣ (3.70) V̂ar = 1 M − k M∑ i−1 (ϕ̂i − ¯̂ϕ)2 (3.71) where M is the number of realization, k is the number of the parameters, ϕ0 is the set of parameter values used to generate the output, ϕ̂ is the set of parameter values estimated from each realization, and the bar denotes the mean value. The estimated bias and variance are presented in Table 3.3. The estimated bias was small and the variance was negligible for the small parameter values. The proposed grey-box identification method thus provides an ideal estimator for the small parameter values. For the large parameter values, the estimated bias and variance became more significant. As Ctot, Ltot and CL increase, K2, K4 and A2 decrease (see equations (3.31),(3.35),(3.36)). Since some of the entries in A(θ) are fixed at 1 (see equation (3.26)), decrease in K2, K4 and A2 makes A(θ) less numerically favorable. This may explain the significant bias and variance for the large parameter values. To avoid this problem, the input and output data can be downsampled before the grey-box identification method is applied. This is illustrated in the next section. 3.2.6 Performance of the proposed grey-box identification method using simulated data To further evaluate the performance of the proposed grey-box identification method, a simulated flow waveform (Fig.3.4) was applied as the input to the model described by equations (3.21) to (3.31). Several sets of parameter values were used as shown in Table 3.4. The input and the model output were downsampled 5 and 10 times. The grey-box identification method was applied to the original and the downsampled input and model output. For the downsampled data, ∆t was increased, and K2, K4 and A2 increased accordingly (see equations (3.31),(3.35),(3.36)). The goal of this simulation study was to decide whether 55 3.2. Methods ϕ0 ˆBias V̂ar Ctot 0.01 0.0005 2.3× 10−10 Ltot 0.01 0.0005 2.8× 10−10 CL 0.02 0.0012 1.4× 10−9 Ctot 0.05 0.0015 2.2× 10−5 Ltot 0.05 0.0020 2.2× 10−5 CL 0.02 0.0007 5.6× 10−6 Ctot 0.10 0.0138 0.0026 Ltot 0.10 0.0730 0.0161 CL 0.02 0.0017 0.0003 Ctot 0.10 0.0092 0.0039 Ltot 0.05 0.0230 0.0012 CL 0.05 0.0014 0.0020 Ctot 0.05 0.0022 0.0003 Ltot 0.10 0.0166 0.0042 CL 0.05 0.0107 0.0219 Table 3.3: The estimated bias and variances: ϕ0 is the set of true parameter values; ˆBias is the estimated bias; V̂ar is the estimated variance; Ctot is the total line capacitance; Ltot is the total line inductance; CL is the load capacitance. 56 3.2. Methods reasonable estimates of K2, K4 and A2, and thereby Ctot, Ltot and CL, can be obtained using the downsampled data. The results are presented in Table 3.4 and 3.5. For the original data, good estimates of K2, K4 and A2, and thereby Ctot, Ltot and CL, were obtained in all cases. For the downsampled data, good estimates of Ctot, Ltot and CL were obtained in all cases. Because of the downsamnpling, ∆t increased 5 or 10 times, and K2, K4 and A2 increased 5 or 10 times accordingly. It is inter- esting to note that, in some cases, the estimated K2, K4 and A2 actually produced an unstable model (Table (3.5)). This was likely caused by the requirement of ∆t imposed by the inequality in (3.11). However, the cor- responding estimates of Ctot, Ltot and CL were still reasonably accurate. Samples of the RAP waveforms estimated using the sets of parameter values identified from the original and the downsampled data are shown in Figure 3.5. For the original data, the error, as evaluated by the differences between the estimated and the true parameter values, increased as the the param- eter values increased (Table 3.4). Similar observation was obtained in the monte-carlo study as described in the previous section. However, for the downsampled data, the error did not seem to be dependent on the parame- ter values (Table 3.4). This is likely due to the increase in K2, K4 and A2, which caused A(θ) to be numerically favorable. On the other hand, as the sampling rate decreased, the error increased (Table 3.4). However, even for the 10-times downsampled data, the error was still acceptable. In this thesis, the sampling rate of the experimental data was similar to the sampling rate of the 10-times downsampled data. 3.2.7 Estimation of the cardiac output from the radial artery pressure waveform Re-estimation of the terminating resistance RT corresponds to the peripheral resistance, and is altered through vasodi- lation and vasoconstriction, which can happens in a very short course of time. Therefore, it is necessary to re-estimate RT when CO is estimated from the identified model during monitoring. To re-estimate RT , the RAP is considered to be entering a Windkessel load (Fig. 3.6). The time constant (τ) of the Windkessel load is given by τ = RTCwind (3.72) 57 3.2. Methods Original sampling rate Downsampled 5 times Downsampled 10 times ∆t = 0.00023 ∆t = 0.0011 ∆t = 0.0023 ϕ0 ϕ̂ ϕ̂− ϕ0 ϕ̂ ϕ̂− ϕ0 ϕ̂ ϕ̂− ϕ0 Ctot 0.05 0.050 0.0000 0.046 -0.0036 0.043 -0.0072 Ltot 0.05 0.050 0.0000 0.053 0.0034 0.058 0.0078 CL 0.02 0.020 0.0000 0.023 0.0034 0.026 0.0064 Ctot 0.05 0.050 -0.0001 0.048 -0.0020 0.0456 -0.0044 Ltot 0.10 0.100 0.0001 0.1039 0.0039 0.1098 0.0098 CL 0.02 0.020 0.0000 0.0219 0.0019 0.0238 0.0038 Ctot 0.10 0.1005 0.0005 0.0967 -0.0033 0.0938 -0.0062 Ltot 0.05 0.0497 -0.0003 0.0506 0.0006 0.0511 0.0011 CL 0.02 0.0201 0.0001 0.0243 0.0043 0.0302 0.0102 Ctot 0.10 0.1025 0.0025 0.0968 -0.0032 0.0928 -0.0072 Ltot 0.10 0.0975 -0.0025 0.1031 0.0031 0.1075 0.0075 CL 0.05 0.0517 0.0017 0.0531 0.0031 0.0566 0.0066 Table 3.4: Physical model parameters estimated from the simulated data: ϕ0 is the set of true parameter values; ϕ̂ is the set of estimated parameter values; Ctot is the total line capacitance; Ltot is the total line inductance; CL is the load capacitance; ∆t is the time step. 58 3.2. Methods Original sampling rate Downsampled 5 times Downsampled 10 times ∆t = 0.00023 ∆t = 0.0011 ∆t = 0.0023 θ0 θ̂ Stablility θ̂ Stablility θ̂ Stablility K4 0.0136 0.0136 0.0732 0.1587 K2 0.0136 0.0136 stable 0.0636 unstable 0.1175 unstable A2 0.0113 0.0113 0.0484 0.0858 K4 0.0136 0.0136 0.0708 0.149 K2 0.0068 0.0068 stable 0.0327 unstable 0.0619 unstable A2 0.0113 0.0113 0.0516 0.0952 K4 0.0068 0.0068 0.0351 0.0724 K2 0.0136 0.0137 stable 0.0671 unstable 0.1330 unstable A2 0.0113 0.0112 0.0466 0.0749 K4 0.0068 0.0066 0.0351 0.0732 K2 0.0068 0.0070 stable 0.0330 unstable 0.0632 unstable A2 0.0045 0.0044 0.0213 0.0400 Table 3.5: Stability of the models identified from the simulated data: K2, K4 and A2 are the unknowns of the physical model; θ0 is the set of the true values of the unknowns at the original sampling rate; θ̂ is the set of the estimated values of the unknowns; ∆t is the time step. 59 3.2. Methods 0 0.05 0.15 0.25 0.35 0.45 −5 0 5 10 15 20 25 30 35 40 45 50 Time (s) F lo w ( cm 3 /s ) Figure 3.4: The simulated flow waveform where Cwind is the Windkessel compliance. τ can be estimated from the pressure waveform using the area method as described by Liu et al. [98]: τ = A12 P1 − P2 (3.73) where P1 is the pressure at the beginning of the diastolic decay, P2 is the end-diastolic pressure, A12 is the area under the pressure waveform between P1 and P2 (Fig. 3.7). During calibration, τ is measured from the RAP waveform, then Cwind is calculated using equation (3.72). When monitoring CO, it is assumed that Cwind is constant. The up- dated τ is measured from the RAP waveform, then RT is re-estimated using equation (3.72). Transfer function from the ascending aortic flow to the radial artery pressure The transfer function (M) from the AAF to the RAP is given by [14]: M(z) = Zch(1 + Γ)z d z2d − Γ (3.74) 60 3.2. Methods 0 0.05 0.15 0.25 0.35 0.45 60 70 80 90 100 110 120 130 Time (s) P re ss u re ( m m H g ) Original D5 D10 Figure 3.5: The estimated radial artery pressure waveforms: original in- dicates the waveform identified from the original data; D5 indicates the waveform identified from the 5-times downsampled data; D10 indicates the waveform identified from the 10-times downsampled data. C wind R T RAP Figure 3.6: The re-estimation of RT by considering the radial artery pressure waveform as entering a Windkessel load: RAP is the radial artery pressure; Cwind is the Windkessel compliance; RT is the terminating resistance. 61 3.2. Methods 0 0.2 0.4 0.6 0.8 1 0 20 40 60 80 100 Time (s) P re ss u re ( m m H g ) P 1 P 2 A 12 Figure 3.7: The determination of the time constant of the Windkessel load using the area method: P1 is the pressure at the beginning of the diastolic decay; P2 is the end-diastolic pressure; A12 is the area under the pressure waveform between P1 and P2. 62 3.2. Methods where z is the z-variable, Γ is the reflection coefficient at the line-load inter- face and d is the discrete time delay between the AAF and the RAP. Since the wave propagation speed for a lossless line is 1/ √ LC, d is calculated as: d = l 1/ √ LC × 1 ∆t = √ (L× l)(C × l) ∆t = p− 1√ K2K4 (3.75) M(z) is given in terms of CL, RL, Zch and RT and d by: M(z) = b1z −d + b2z −d−1 1 + a1z−1 + a −2d 2 + a −2d−1 3 (3.76) where b1 = 2Zch(1 + 2FsCLRL) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (3.77) b2 = 2Zch(1− 2FsCLRL) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (3.78) a1 = Zch RT + 1− 2FsCL(RL + Zch + RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (3.79) a2 = Zch RT − 1− 2FsCL(RL − Zch − RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (3.80) a3 = Zch RT − 1 + 2FsCL(RL − Zch − RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (3.81) where Fs is 1/∆t. The derivation of equation (3.76) is given in Appendix C. AAF is estimated from RAP using the inverse of M(z). Calculation of cardiac output from the model-estimated ascending aortic flow waveform CO can be calculated as the area integral (AL) of the model-estimated AAF waveform divided by its duration (D). During calibration, CO is also estimated using echocardiography. A correction factor (β) is derived from AL and D of the model-estimated AAF waveform and CO estimated using echocardiography (COecho): β = COecho AL/D (3.82) 63 3.3. Experiments During monitoring, CO estimated from the RAP (CORAP) is calculated as: CORAP = β AL D (3.83) Summary of the CO estimation technique The steps of the CO estimation technique are shown in Figure 3.8. 3.3 Experiments 3.3.1 Data acquisition This study was approved by the University of British Columbia Children’s and Women’s Research Ethics Board. Five children who underwent car- diac surgeries participated. Left ventricular outflow tract (LVOT) diameter was obtained using two-dimensional echocardiography. Flow velocity at the LVOT was obtained using pulsed-Doppler ultrasound (Fig. 3.9). Pressure was obtained at the radial artery from an indwelling catheter (Fig. 3.10). Velocity and pressure data were obtained simultaneously. The echocardio- graphic data were obtained using a GE Vivid 7 Pro ultrasound machine (GE Vingmed, Horten, Norway). The pressure data were obtained using a Philips IntelliVue MP70 Patient Monitor (Philips Healthcare, Andover, MA, the United States). Data were collected once the patients were deemed stable after the surgeries (0 h). In 4 patients, data were also collected 3 to 5 hours after (5 h). In the remaining patient, data were collected at the beginning (0 h) and 24 hours after (24 h). 3.3.2 Data processing The envelope of the Doppler spectrum of the LVOT flow velocity was traced automatically (Fig. 3.9). The LVOT cross-sectional area was calculated from the LVOT diameter. The LVOT flow waveform was obtained by mul- tiplying the LVOT flow velocity with the LVOT cross-sectional area. The LVOT flow waveform was used as a surrogate of the AAF waveform. The pressure waveform showed a cyclic fluctuation with a period of about 2 to 3 second (Fig. 3.10). It was assumed that this fluctuation was due to respiration. For the calibration and each subsequent CO measurement, five to ten consecutive beats were averaged. The number of beats averaged was chosen to include one complete respiratory cycle (Fig. 3.10). 64 3.3. Experiments R T is estimated from I 1 and V p-1 using equation (3.32) L tot ,C tot , R L , C L are identified from I 1 and V p-1 using greybox identification C wind is estimated from V p-1 and R T using equations (3.54) and (3.55) β is derived from the model-estimated AAF and the measured echocardiographic cardiac output using equation (3.65) I 1 is the measured AAF and V p-1 is the measured RAP C al ib ra ti o n R T is re-estimated from the measured RAP and C wind using equations (3.54) and (3.55) M(z) is defined from R T , L tot ,C tot , R L , C L using equations (3.58) to (3.64) Model-estimated AAF is obtained from the measured RAP using the inverse of M(z) CO is estimated from the model-estimated AAF and β using equation (3.66) E st im at io n Figure 3.8: The cardiac output estimation technique: I1 is the current at the beginning of the line; AAF is the ascending aortic flow; Vp−1 is the voltage at the end of the line; RAP is the radial artery pressure; RT is the terminating resistance; Ltot is the total line inductance; Ctot is the total line capacitance; RL is the load resistance; CL is the load capacitance; Cwind is the Windkessel compliance as seen by the radial artery pressure; β is a correction factor relating the cardiac output estimated by echocardiography and that estimated from the radial artery pressure; M(z) is the transfer function from the ascending aortic flow to the radial artery pressure; CO is the cardiac output. 65 3.3. Experiments −0.2 0.2 0.6 1.0 1.4 0 0.5 1 1.5 2 2.5 Time (s) V el o ci ty ( m /s ) Figure 3.9: The Doppler spectrum of the LVOT flow velocity: the envelope of the spectrum was traced automatically. For the pressure waveform, the beginning and end of each beat were marked by the end-diastolic pressures. The beginning and end of the chosen pressure beat were equalized by removing the linear component defined by these two points. For the flow waveform, only the positive flow during the ventricular ejection was considered. Zero flow was assumed after the ventricular ejection (Fig. 3.11). The sampling rate of the flow waveform was in the range of 455 – 500 Hz (∆t was 0.0020 – 0.0022s). The sampling rate of the pressure waveform was 125 Hz (∆t was 0.0080s). The pressure waveform was upsampled to have the same sampling rate as the flow waveform. 3.3.3 Cardiac output determination CO was estimated from the RAP waveform using the proposed technique. To evaluate the accuracy of the proposed technique, CO was also estimated using echocardiography. Using echocardiography, each measurement of CO was calculated as the velocity-time-integral of the LVOT flow velocity waveform times the LVOT area divided by the duration of the waveform. Each determination of CO was the average of one to three measurements of CO. 66 3.3. Experiments 1 2 3 4 5 6 7 8 9 10 30 40 50 60 70 80 90 Time (s) P re ss u re ( m m H g ) 1 respiratory cycle Figure 3.10: The radial artery pressure waveform: one respiratory cycle is indicated. 0 0.1 0.2 0.3 0.4 −10 0 10 20 30 40 50 60 Time(s) F lo w ( cm 3 /s ) Figure 3.11: The averaged ascending aortic flow waveform: zero flow was assumed after the ventricular ejection. 67 3.4. Results Patient Age (mo) BSA (m2) 1 1 0.24 2 2 0.27 3 5 0.33 4 6 0.34 5 9 0.30 Table 3.6: Demographics of patients: BSA is the body surface area. 3.4 Results The demographics of the patients are listed in (Table 3.6). The system order was determined during subspace identification. Singu- lar value decomposition was performed on a matrix that has the same column space as the extended observability matrix. The number of dominant singu- lar values corresponds to the system order. In a previous modelling study on simulated data by Kind et al. [85], it was found that for noise-free simu- lated data, the gap between the dominant and non-dominant singular values was easily distinguishable. However, as noise was added, this gap narrowed. For the experimental data in this thesis, the singular value decreased rather smoothly, and the gap between the dominant and non-dominant singular values was not obvious (Fig 3.12). To determine the system order, it was assumed that the physical model should have at least two elements, which corresponds to a system order of 6. A relatively significant drop of the log of the singular value could be identified between the system order of 6 to 12. It was observed that the drop usually occurred from 1 to 0. There- fore, the system order was chosen as the maximum system order with the corresponding log of singular value greater than 0.1. The identified physical parameters are listed in Table 3.7. In all five cases, the error criterion and the model parameters converged rapidly during the iterations of bilinear regression. Samples of the convergence graphs are shown in Figure 3.13. The CO determinations and the corresponding time of measurements are listed in Table 3.8. To evaluate the accuracy of our proposed technique, echocardiography was used as a reference method. CORAP is the CO es- timated from the RAP waveform using equation (3.83). COecho is the CO estimated using echocardiograhy. The bias and limits of agreement are pre- sented in terms of % error (Table 3.9), CO (Table 3.10), and cardiac index (CI) (Table 3.11) for the data collected during different time intervals. CI is 68 3.4. Results 2 4 6 8 10 12 −2 0 2 4 6 8 10 12 L o g o f si n g u la r v al u e System order 0.91 -0.24 Figure 3.12: Log of singular value versus system order Patients Ctot (cm3/mmHg) Ltot (mmHgs2/cm3) CL (cm3/mmHg) RL (mmHgs/cm3) 1 0.0887 0.0559 0.0185 0.9904 2 0.0350 0.0765 0.0083 1.6793 3 0.0606 0.0597 0.0450 1.2892 4 0.0854 0.0530 0.0212 0.9490 5 0.0325 0.1356 0.0219 3.0741 Table 3.7: Identified physical model parameters: Ctot is the total line ca- pacitance; Ltot is the total line inductance; CL is the load capacitance; RL is the load resistance. 69 3.5. Discussion Patient Time of measurements COecho (L/min) CORAP (L/min) % Difference 1 0h 0.924 0.877 5.3 5h 0.983 0.937 4.8 2 0h 0.253 0.247 2.4 3h 0.347 0.299 14.7 3 0h 0.776 0.754 2.8 4h 0.994 0.979 1.5 4 0h 0.828 0.744 10.7 5h 0.826 0.885 -6.9 5 0h 0.638 0.633 0.8 24h 0.791 0.559 34.4 Table 3.8: Cardiac output determinations: COecho is the cardiac output de- termination obtained using echocardiography; CORAP is the cardiac output determination obtained from the radial artery pressure waveform. Time of measurements Bias (%) Limits of agreement (%) 0 hrs 4.4 -3.2 – 12.0 0, 5 hrs 4.0 -8.0 – 16.1 0, 5, 24 hrs 7.1 -14.9 – 29.0 Table 3.9: Bias and limits of agreement in terms of % error CO normalized to the body surface area (BSA), and it is also often measured during critical care. The Bland-Altman plot is shown in Figure 3.14. A sample of the model- estimated AAF waveform and the corresponding measured AAF waveform are shown in Figure 3.15. 3.5 Discussion In this study, the reference method was echocardiography. Calibration was performed only at the beginning of the data collection. For a 5-hour calibration-free period, the limits of agreement was within ±17%. To put this into perspective, studies of other pulse contour techniques were consid- ered. There are two relatively well-established techniques based on the pulse contour analysis [67]: PiCCO [47], [54], [64] and PulseCO [63], [84]. Many 70 3.5. Discussion 0 500 1000 1500 2000 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 0 500 1000 1500 2000 0.1 0.105 0.11 0.115 0 500 1000 1500 2000 0.086 0.088 0.09 0.092 0.094 0.096 0.098 0.10 0 500 1000 1500 2000 0.04 0.05 0.06 0.07 0.08 0.09 0.1 e K 4 K 2 A 2 Iteration Iteration Iteration Iteration Figure 3.13: Samples of convergence graphs: e is the error criterion; K2, K4 and A2 are the unknowns of the physical model. Time of measurements Bias (L/min) Limits of agreement (L/min) 0 hrs 0.033 -0.033 – 0.099 0, 5 hrs 0.024 -0.055 – 0.103 0, 5, 24 hrs 0.045 -0.104 – 0.194 Table 3.10: Bias and limits of agreement in terms of cardiac output Time measurements Bias (L/min/m2) Limits of agreement (L/min/m2) 0 hrs 0.110 -0.098 – 0.319 0, 5 hrs 0.088 -0.169 – 0.344 0, 5, 24 hrs 0.156 -0.333 – 0.646 Table 3.11: Bias and limits of agreement in terms of cardiac index 71 3.5. Discussion 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 −50 −40 −30 −20 −10 0 10 20 30 40 50 Averaged CO (L/min) % D if fe re n ce mean+1.96SD mean mean-1.96SD Figure 3.14: Bland-Altman plot of the cardiac output determinations ob- tained using the echocardiography technique and from the radial artery pres- sure waveform; CO is the cardiac output. 72 3.5. Discussion 0 0.5 1 1.5 2 2.5 −10 0 10 20 30 40 50 60 F lo w ( cm 3 /s ) Time (s) Measured flow Estimated flow Figure 3.15: Measured and estimated ascending aortic flow waveforms studies have been conducted to evaluate the accuracy of these techniques on different patient populations under different conditions. To make a reason- able comparison, studies conducted on similar patients and/or under simlar conditions were considered. In a study of PiCCO [47], CO determinations were obtained from patients who had underwent pediatric cardiac surgeries. Calibration was performed every 4 hours. The limits of agreement of CI were found to be approximately ±0.75 L/min/m2. In a study of PulseCO [84], CO determinations were obtained from children undergoing hemodynamic evaluation in the cardiac catheterization laboratory. Data were collected in a 4-hour calibration-free period. The limits of agreement of CI were found to be approximately -0.09–0.47 L/min/m2. The accuracy of our technique, in terms of limits of agreement, is comparable to the accuracies of PiCCO and PulseCO. It has been suggested that the measurement of physiological variables often lacks precision and has an error of ±10–20% [36]. For example, it has been shown that thermodilution has an inherent error (variability) of 13% (average of three measurements) and 22% (single measurements) [152]. Therefore, to evaluate a new CO estimation technique with thermodilution as the reference technique, the limits of agreement are expected to exceed 73 3.5. Discussion ±20 %. Critchley and Critchley [36] then proposed that limits of agreement within ±30% are acceptable, by assuming both of the test and reference techniques have an error of 20%. In this study, the limits of agreement including all CO determination (0h, 5h, 24h) were indeed within ± 30%. However, it was noted that, the % difference was much larger for the CO determination obtained at 24 h (Fig. 3.14). The limits of agreement excluding the last measurement (0h, 5h) were only ±17%, and can be considered as quite ideal. It suggests that a re-calibration after 5 hours may improve the performance of the technique. A major limitation of this study is that only two CO determinations were obtained from each patient. Therefore, the performance of the proposed technique over time cannot be evaluated . It would be necessary to determine how long the technique would stay accurate without re-calibration. Both PiCCO and PulseCO require calibration using transpulmonary lithium-dilution or thermodilution. Transpulmonary lithium-dilution re- quires the placement of both central venous and arterial catheters. A ma- jor advantage of our technique is that the calibration is completely non- invasive. Since the ultrasound machine has become more compact and portable, echocardiographic examinations can be conducted in ICU with relative ease. The proposed algorithm can be implemented as a module with data transfer interfaces to the ultrasound machine and the bedside monitor. Recently, it has been suggested that central hemodynamics can be es- timated from multiple peripheral pressures [60], [107], [155]. Cardiac out- put estimation is thus possible using multiple peripheral pressures that can be obtained minimally-invasively or non-invasively. However, validations of these techniques have yet to be performed on human data. 74 Chapter 4 Measuring Aortic Pressure Waveform using Ultrasound B-Mode Imaging1,2 Aortic pressure is closely related to the pathogenesis of cardiovascular dis- eases. A non-invasive method for the measurement of the aortic pressure waveform would be a very valuable tool for the diagnosis of cardiovascular diseases at their early stage of development. In this chapter, a method for estimating the aortic pressure waveform using B-mode imaging is proposed. A semi-automatic algorithm is developed to measure the aortic distension waveform from the B-mode images of the aortic arch acquired at a supraster- nal view. A linear calibration scheme and an exponential calibration scheme are used to estimate the aortic pressure waveform from the aortic distension waveform. 4.1 Introduction The aortic pressure waveform is different in both amplitude and morphol- ogy from the peripheral artery pressure waveform [119]. To assess the aor- tic pressure waveform non-invasively, the carotid artery pressure waveform obtained using applanation tonometry can be used as a surrogate to the aortic pressure waveform [82]. However, the use of applanation tonometry on the carotid artery requires substantial expertise. Recently, there has been increasing interest in estimating the pressure waveform from the dis- tension waveform, which is the change in the arterial diameter over time. Various methods have been proposed to obtain the distension waveform 1 A version of this chapter is under preparation for submission. 2 M. Leung, G. A. Dumont, G. G. S. Sandor, M. T. Potts and J. E. Potts, “A novel method to estimate the aortic pressure waveform using B-mode ultrasound images acquired from a suprasternal view,” 30th IEEE Engineering in Medicine and Biology Conference, Vancouver, 2008 75 4.1. Introduction using M-mode imaging [151], B-mode imaging [55,56,151], echo-tracking [86,151,153,161,166], and more recently, tissue Doppler imaging [41,61]. In medical ultrasound imaging, pulses of sound are sent into the body. The sound wave propagates through the tissue and is scattered by small changes in density. The backscattered signal is received by the transducer. Its strength is determined by the deviations within the tissue and the com- position of the tissue, and its transit time is determined by the depth of the tissue. Therefore, the backscattered signal received along a scan line provides the information of the tissue layers along that line. In B-mode imaging, “B” stands for brightness. The scan line sweeps through a plane. The amplitude of the backscattered signal is displayed as a grey value. The backscattered signals received along the scan line as it sweeps through the plane are displayed side-by side, and a two-dimensional image is thus con- structed. In M-mode imaging, “M”stands for motion. Sound waves are emitted successively along a scan line. The amplitude of the backscattered signal is also displayed as a grey value. The backscattered signals received at successive times are displayed side-by-side, and thereby generates an image of tissue motion at the direction of the scan line. The principle of echo- tracking is quite different from that of B-mode and M-mode imaging. In echo-tracking, pulses of sound are emitted successively, and the returned signal from a particular echo structure, e.g. arterial wall, is identified by its amplitude or phase. The transit time of this returned signal provides infor- mation on the depth of the structure, and the movement of the structure can be obtained from the change in the transit time. Tissue Doppler imaging is based on the same principle as the conventional Doppler echocardiography, which relies on the detection of the shift in frequency of the returned signal from moving objects. In tissue Doppler imaging, the Doppler principle is used to quantify tissue motion. Of all these methods, echo-tracking has the advantage of a high spatial resolution of about 1 µm [66]. In contrast, B-mode imaging can only achieve a moderate spatial resolution of about 0.3 mm [151]. However, B-mode imaging also provides additional spatial information along the artery. It has been suggested that,with the use of proper mathematical modelling, this additional spatial information can be used to obtain the distension waveform [150]. It has been shown that distension waveform at the carotid artery obtained using B-mode imaging with modelling has similar reproducibility as that obtained using echo-tracking [151]. In this chapter, an algorithm for obtaining the aortic distension waveform at the ascending aorta from the B-mode longitudinal images of the aortic arch is developed. In addition to the modelling of the spatial structure 76 4.1. Introduction of the aorta, this new algorithm also makes use of the continuity of the aortic wall motion to obtain a more accurate distension waveform. This new algorithm consists of three steps. First, an adaptive block matching scheme is used to identify the approximate locations of the edges of the aortic arch. Second, the approximate edges are refined according to the gradient of the intensity of the image. Second-degree polynomials are then fitted to the refined edges. Third, the refined edges are further fine-tuned by imposing smoothness constraints on the curvatures and vertexes over time. The distension waveform is the distance between the edges over a cardiac cycle. To estimate the pressure waveform from the distension waveform, a re- lationship between the pressure and the diameter has to be assumed. Even though the arterial wall is viscoelastic and exhibits nonlinear elasticity, it has been shown that the pressure waveform and the distension waveform at the ascending aorta of human have gross similarity [57]. Therefore, a linear relationship may be assumed between the pressure and the diameter. The distension waveform may also be estimated from the pressure waveform through an exponential calibration scheme. The relationship between the pressure and the cross-sectional area of the artery can be described by an exponential function [110]. The cross-sectional area of the artery can be found from its diameter. The exponential function involves a wall rigidity coefficient, which can be found through an iterative scheme [110]. In this chapter, an aortic pressure waveform measurement method is developed. In this method, the aortic distension waveform is obtained from B-mode images of the aortic arch using the semi-automatic edge detection algorithm developed in this chapter, then the aortic pressure waveform is derived from the aortic distension waveform using an exponential calibration scheme proposed by previous researchers [110]. To investigate the clinical usefulness of the aortic pressure waveform mea- surement method developed in this chapter, a clinical study was performed in which this new method was applied to 9 healthy children. Aortic pressure waveforms were estimated from the aortic distension waveform through both linear and exponential calibration schemes. Carotid artery applanation to- mometry was used as a reference method. Carotid artery pressure waveform was obtained using applanation tonometry and used as a surrogate of the aortic pressure waveform. The tonometric and B-mode waveforms were com- pared in the time domain through the root-mean-squared differences, and in the frequency domain through the moduli and phases of the harmonics. The aortic pulse pressure (PP) and two arterial stiffness indices, augmentation index (AI) and total arterial compliance (TAC), were calculated from the 77 4.2. Methods Figure 4.1: Left: several kernel blocks are selected along the edges on the first image. Right: a best-matching block is found inside the search region on the next image. waveforms. 4.2 Methods 4.2.1 Adaptive block matching In this semi-automatic edge detection algorithm, the aortic distension wave- form is measured from a series of B-mode images of the aortic arch collected over one cardiac cycle. On the first image, several points along the edges are selected manually (Fig. 4.1). These selected points give the approximate locations of the edges. Using block matching, these approximate locations are tracked on the subsequent images. Kernel blocks of size k × k are defined with the selected points as the centers of the blocks. For each kernel block, a search region of size s × s is defined with the kernel block at the center (Fig. 4.1). On the subse- quent images, the best-matching block to the kernel block is found within the search region. The best-matching block has the maximum correlation with the kernel block. The correlation (cr) between the kernel block and a candidate block is given by cr = ∑k x=1 ∑k y=1K(x, y)× C(x, y)√∑k x=1 ∑k y=1[K(x, y)] 2 ×∑kx=1∑ky=1[C(x, y)]2 (4.1) where K is the intensity of the kernel block and C is the intensity of the candidate block, and x and y are the spatial indices. 78 4.2. Methods Frame 1 Frame 20 Figure 4.2: Substantial change of the curvature of the aortic wall: the part of the edge inside the original kernel block on the first image (yellow) becomes very different on the twentieth image (green); the best-matching block (red) to the original kernel block does not contain the same part of edge. As the curvature of the aortic wall changes substantially through out the cardiac cycle, the edge contained by the kernel block becomes much different. The candidate block that has the maximum correlation with the kernel block does not necessarily contain the same edge (Fig. 4.2). To track the edge during pulsation, an adaptive scheme is used. When the maximum correlation does not reach a certain threshold, the best-matching block of the previous image is used as the kernel block instead. Since the aortic wall motion is gradual, using the previous best-matching block as the kernel block can track the changes in the wall curvature (Fig. 4.3). However, if the previous best-matching block is used as the kernel block on all images, the information in the original kernel block may be completely lost and the best-matching block may eventually become too deviated from the original kernel block. Therefore, the previous best-matching block is used only when the maximum correlation does not reach a certain threshold γ, chosen to be 0.95. In the experiments described in section 4.3, the image size was approx- imately 400 by 600 pixels. The actual size of the image was approximately 79 4.2. Methods Frame 1 4 7 10 13 16 19 22 25 28 31 34 37 40 43 Figure 4.3: The adaptive block matching scheme: the set of best-matching blocks to the original kernel block (shown in figure 4.2) are found by the adaptive block matching scheme; the cyan highlighted blocks are found by comparing with the best-matching block on the previous image; the non- highlighted blocks are found by comparing with the original kernel block. 10 by 15 cm. k was chosen to be 4 mm (approximately 17 pixels) so that the edge contained by the kernel block had trackable movement. s was chosen to be 10 mm (approximately 43 pixels) so that the edge contained by the kernel block was always with in the search region. 4.2.2 Gradient-based edge refining A preliminary edge is defined by joining the centers of the best-matching blocks using a cubic interpolating spline [131]. The more accurate loca- tions of the edge points are then searched for along the normal direction of the preliminary edge (Fig.4.4). The intensity function along each normal is smoothed using a cubic smoothing spline [131]. The location of the maxi- mum gradient of the smoothed intensity function is the refined location of the edge point (Fig.4.4). The edge refining scheme gives a set of refined edge points. For the outer edge, some of these edge points are located at the “disconnections”. These “disconnections” are due to the presence of the innominate artery, innominate vein, the left common carotid and the left subclavian arteries. These edge points do not represent the actual edge and are excluded by comparing with the preliminary edge. The discrepancy is defined as the norm between each preliminary edge point and the corresponding refined edge point. If the discrepancy is more than β, the edge point is excluded. 80 4.2. Methods 20 40 60 80 100 120 140 160 180 200 In te n si ty Figure 4.4: Left: the more accurate locations of the edge points are searched for along the normals of the preliminary edge. Right: a smoothing spline (green dash line) is fitted to the intensity function (blue solid line) of the normal highlighted in red on the image on the left; the vertical red line indicates the refined location of the edge. In the experiments described in section 4.3, β was chosen to be 12 pixels for the outer edge and 8 pixels for the inner edge. Similar to the choice of k and s, β was chosen based on the resolution of the image and the normal anatomical size of the aorta and the adjacent vessels. 4.2.3 Weighted total least-squares fitting A second degree polynomial given by: y = ax2 + bx+ c (4.2) is fitted to the refined edge points. Since the refined edge points are found by searching along the normal direction of the preliminary edge, the error of the fit is quantified by the perpendicular distances between the points and the fit. Total least-squares fitting is used to minimize this error [163]. During total least-squares fitting, weights are assigned to the refined edge points. Two weight assigning schemes are used. The first scheme is based on the assumption that a refined edge point is located either on the actual edge or some unrelated bright speckles. If a refined edge point is indeed located on the actual edge, it should be adjacent to a few other edge points. For each refined edge point, the weight assigned is directly proportional to the number of other refined edge points inside its surrounding region. If the number of points is less than or equal to θ, zero 81 4.2. Methods Figure 4.5: Left: some of the refined edge points are located at the edges of the left common carotid and the left subclavian arteries. Right: the total least-squares fit obtained without the second weight assigning scheme (red solid line) has a higher curvature than the total least-squares fit obtained with the second weight assigning scheme (cyan dash line). weight is assigned. If the number of points is more than θ, the assigned weight is the number of points times 10. This weight-assigning scheme is applied to both outer and inner edges. In the experiments described in section 4.3, θ was chosen to be 5 for the outer edge and 3 for the inner edge in this thesis. θ was chosen based on the resolution of the image and the normal anatomical size of the aorta and the adjacent vessels. In some images, some of the refined edge points are located at the edges of the left common carotid or left subclavian arteries (Fig.4.5). These refined edge points cause the total least-squares fit to have a higher curvature. The influence of these refined edge points are reduced by assigning to them a smaller weight. The refined edge points are ordered clockwise. For each refined edge point, its distances with the next three refined edge points are calculated. Its weight is inversely proportional to the sum of these distances. When the sum of distances is over ϕ, a very small weight that is practically zero is assigned. This weight assigning scheme is applied only to the outer edge. In the experiments described in section 4.3, ϕ was chosen to be 50. ϕ was chosen based on the resolution of the image and the normal anatomical size of the aorta and the adjacent vessels. 82 4.2. Methods 4.2.4 Continuity of the aortic wall motion During systole, the aorta expands in order to accommodate the sudden in- crease in blood volume, and then slowly returns to its previous state through- out diastole. The aortic wall motion is gradual. It is therefore assumed that the changes in the curvature and the vertex of the total least-squares fit over time are smooth. The curvature (cv), the x-coordinate of the vertex (vx) and the y-coordinate of the vertex (vy) of a second degree polynomial are given by: cv = a (4.3) vx = −b 2a (4.4) vy = −(b2 − 4ac) 4a (4.5) The functions of cv, vx and vy over time are smoothed by a cubic smoothing spline [131]. The finalized edges are then obtained. 4.2.5 Choosing a cross-section along the arch The cross-section is chosen at the ascending aorta proximal to the innomi- nate artery, and as perpendicular to both the outer and inner edges as pos- sible. Three cross-sections are chosen. The distension waveforms obtained along each cross-section are averaged to give the final distension waveform. 4.2.6 Estimating the pressure waveform from the distension waveform Mean arterial pressure and diastolic pressure remain relatively constant through out the arterial tree [89]. It is thus assumed that the mean and diastolic pressures obtained at the brachial artery are the same as the mean and diastolic pressures at the aorta. Systolic pressure (SP) and diastolic pressure (DP) are measured at the upper arm using a cuff sphygmomanometer. Mean pressure (MP) is calcu- lated from the traditional formula [21]: MP = 1 3 SP + 2 3 DP (4.6) 83 4.3. Experiments Linear calibration It is assumed that the pressure and diameter are linearly related. DP is assigned to the minimum of the distension waveform, and MP is assigned to the area integral divided by the duration of the distension waveform. Exponential calibration It is assumed that the pressure and the cross-sectional area of the artery are related by an exponential function [110]: P (t) = P0 exp α( A(t) A0 −1) (4.7) where P (t) is the pressure function, P0 is the reference pressure, A(t) is the cross-sectional area function, A0 is the cross-sectional area at P0, and α is the wall rigidity coefficient. The cross-section area is related to the diameter (D) by A = (D/2)2π. DP is assigned to P0, then the minimum of the distension waveform corresponds to A0. α is found through an iterative search. With an initial guess of α, the mean aortic pressure (P̄a) can be calculated from (4.7). If the estimated P̄a is different from MP, α is recalculated as: αnew = MP P̄a α (4.8) and then P̄a is re-estimated. The iterative search stops when P̄a and MP differ by less than 0.1 mmHg. 4.3 Experiments 4.3.1 Data acquisition and processing 9 healthy children participated in this study. B-mode images of the aortic arch were obtained from a suprasternal view using a 2.2-5.0 MHz phased-array transducer. The transmittal frequency was 2.4 MHz. The frame rate ranged from 49 to 63 frames per second. The carotid artery pressure waveform was obtained using applanation tonometry. The ascending aortic blood velocity was obtained simultaneously using Doppler ultrasound with the transducer positioned at the suprasternal notch. 84 4.3. Experiments The aortic annulus was measured using two-dimensional echocardiogra- phy. Brachial artery systolic and diastolic blood pressures were measured using a cuff sphygmomanometer. All the applanation-tonometric and echocardiographic data were ac- quired using a GE Vivid 7 Pro ultrasound machine (GE Vingmed, Horten, Norway) equipped with an applanation tonometer. Aortic flow was calculated as the envelope of the Doppler spectrum of the aortic blood velocity times the lumen area. The lumen area was assumed to be circular and was calculated from the aortic annulus. The mean aortic flow was used to calculate the peripheral resistance. 4.3.2 Calibration of the carotid artery pressure waveform obtained using applanation tonometry To calibrate the carotid artery pressure waveform obtained using applana- tion tonometry, DP is assigned to the minimum of the waveform and MP is assigned to the area integral divided by the duration of the waveform. 4.3.3 Comparison of waveforms The waveforms were compared in the time domain using the root-mean- squared difference (RMSD). RMSD is given by √∑ i (Ptono(i)− Pdw(i))2 n (4.9) where Ptono is the pressure waveform obtained using applanation tonome- try, Pdw is the pressure waveform estimated from the distension waveform using the linear or the exponential calibration scheme, and n is the number of data points.Pdw was first upsampled to have the same data points as Ptono using a cubic interpolating spline [131]. The waveforms were also compared in the frequency domain using the moduli and the phases of the first three harmonics of the waveforms. 4.3.4 Calculations of pulse pressure, augmentation index and total arterial compliance PP was calculated as the peak of the pressure waveform (Pmax) minus DP. To calculate AI and TAC, Pdw was first upsampled to have the same data points as Ptono using a cubic interpolating spline [131]. 85 4.4. Results Age (yrs) 11± 2 Gender (M:F) 3 : 6 BSA (m2) 1.2± 0.2 SP (mmHg) 108± 7 DP (mmHg) 64± 7 AO (cm) 1.65 ± 0.14 Table 4.1: Subjects’ demographics: BSA is the body surface area; SP is the systolic pressure at the upper arm; DP is the diastolic pressure at the upper arm; AO is the aortic annulus. AI was calculated as the ratio of the difference between Pmax and the pressure at the inflection point (Pinf) to PP: AI = Pmax − Pinf PP (4.10) If Pinf occurs after Pmax, AI is negative [116]. In a study by Segers et al. [142], the inflection point was identified as the middle of the two zero-crossings of the second derivative of the pressure waveform. The zero- crossing of the second derivative corresponds to the local extremum of the first derivative. In this thesis, the inflection point was identified as the mid- dle of a local minimum and a local maximum of the first derivative. This method is equivalent to the method used by Segers et al. TAC was calculated using the area method [98]: TAC = A12 RP (P1 − P2) (4.11) where P1 is the pressure at the very beginning of the diastolic decay, P2 is the end-diastolic pressure, A12 is the area under the pressure curve between P1 and P2, and RP is the total peripheral resistance. RP was estimated as the mean pressure divided by the mean flow. 4.4 Results Data are expressed as mean±SD. A value of p < 0.01 was considered as statistically significant. The demographics of the subjects are listed in Table 4.1. Plin is the aortic pressure waveform obtained from the aortic distension waveform through the linear calibration scheme. Pexp is the aortic pres- 86 4.4. Results sure waveform obtained from the aortic distension waveform through the exponential calibration scheme. The wall rigidity coefficient was 1.2± 0.5 4.4.1 Waveform comparisons In all 9 subjects, Plin and Pexp were compared with Ptono, which was used as the reference. RMSD between Ptono and Plin was 3.9 ± 1.5 mmHg. RMSD between Ptono and Pexp was 4.1 ± 1.5 mmHg. The difference in RMSDs obtained by the two calibration schemes was statistically insignificant. Therefore, in terms of RMSD, the two calibration schemes performed equally. The moduli and phases of the first three harmonics of Ptono, Plin and Pexp are listed in Table 4.2. For the first and second harmonics, the moduli of Ptono were only moderately correlated with that of Plin and Pexp. How- ever, in terms of the absolute pressure, the differences were quite small. The differences were also statistically insignificantly. For the third harmonic, the modulus of Ptono was highly correlated with that of Plin and Pexp. Again, in terms of the absolute pressure, the differences were very small. For the modulus of the third harmonic, the difference between Ptono and Pexp was statistically insignificant, but the difference between Ptono and Plin was statistically significant. The differences in phase were statistically insignifi- cant, except for the first harmonic between Ptono and Plin. Overall, in terms of the differences in the frequency components, the exponential calibration scheme performed better than the linear calibration scheme. 4.4.2 Comparisons of pulse pressures, augmentation indices and total arterial compliances PP, AI and TAC calculated from Ptono, Plin and Pexp in all 9 subjects are listed in Table 4.3. PP estimated from Ptono was highly correlated with that estimated from Plin and Pexp. The differences were small and statistically insignificant. AI estimated from Ptono was moderately correlated with that estimated from Plin and Pexp. AI estimated from Plin was significantly different from that estimated from Ptono. In contrast, the difference in AIs estimated from Ptono and Pexp was small and statistically insignificant. TAC estimated from Ptono was highly correlated with that estimated from Plin and Pexp. The differences were small and statistically insignif- icant. However, the difference obtained from the exponential calibration 87 4.4. Results Harmonics Ptono Plin Pexp ∆tono−lin ∆tono−exp Modulus (mmHg) 1st 5.6± 1.3 6.2 ± 1.4 6.8± 1.5 −0.6± 1.3 −1.2 ± 1.4 r = 0.49 r = 0.52 p =NS p =NS 2nd 3.0± 0.7 2.5 ± 0.6 2.8± 0.7 0.5± 0.6 0.2± 0.6 r = 0.60 r = 0.58 p =NS p =NS 3rd 1.7± 0.8 1.3 ± 0.9 1.4± 1.1 0.4± 0.4 0.3± 0.5 r = 0.91 r = 0.92 p < 0.009 p =NS Phase 1st −1.8± 0.2 −2.0± 0.1 −1.9± 0.1 0.1± 0.1 0.1± 0.1 r = 0.73 r = 0.76 p < 0.005 p =NS 2nd −2.2± 0.3 −2.4± 0.4 −2.5± 0.5 0.1± 0.4 0.2± 0.4 r = 0.48 r = 0.53 p =NS p =NS 3rd −1.6± 2.5 −2.1± 1.9 −0.3± 3.0 0.6± 1.7 −1.3 ± 2.6 r = 0.73 r = 0.56 p=NS p =NS Table 4.2: Comparisons of the waveforms obtained in 9 healthy children using applanation tonometry and B-mode imaging in the frequency domain: Ptono is the carotid artery pressure waveform obtained using applanation tonometry; Plin is the aortic pressure waveform obtained using B-mode imaging with linear calibration; Pexp is the aortic pressure waveform ob- tained using B-mode imaging with exponential calibration; ∆ indicates the difference; r is the correlation between Ptono and Plin, or Ptono and Pexp. 88 4.5. Discussion scheme was much smaller than that obtained from the linear calibration scheme. 4.5 Discussion The pressure-diameter relationship of the arterial wall is complicated by two aspects: nonlinear elasticity, which is often studied as static elasticity, and viscoelasticity, which is often studied as dynamic elasticity. The nonlinear elasticity of the arterial wall is due to its complex com- position. The predominant elastic materials of the arterial wall are elastin and collagen [119]. As the distending pressure inside an artery increases, the circumferential stress transfers progressively from the elastin on to the relatively inextensible collagen, and results in a nonlinear stress-strain, or increasing elastic modulus-pressure, relationship [12]. In the in vitro study of static elasticity of human arteries by [95], the circumferential Young’s modulus of the upper thoracic aorta increased slightly and slowly up to 100 mmHg, and then it began to increase more rapidly. It implies that the pressure-diameter relationship is relatively linear up to 100 mmHg. The elastic behavior caused by step pressure is rather different from that caused by the rapid pressure change during each heartbeat [13]. This is because of the viscoelasticity of the arterial wall. When a viscoelastic material is extended rapidly, it is stiffer than when it is extended slowly [119]. It is said to have a retarded elastic response to an applied stress as it will slowly extend to its final length. When the stress is removed, the reverse phenomenon occurs. The pressure-diameter relationship in human has been studied in vivo [72]. The diameter was found to have a phase lag of -6.7◦ to the pressure at the fundamental frequency of 1.2 Hz. This corresponds to a time delay of 0.016s. For a normal heart rate of 70 bpm, a delay of 0.016s is about 2 % of the cardiac cycle. Combining the effects of nonlinear elasticity and viscoelasticity, the dis- tension waveform should have similar morphology as the pressure waveform during normal pulsation. In this study, it was found that the pressure wave- form and the distension waveform indeed had similar morphology as indi- cated by the low RMSD between Ptono and Plin. The morphologies of Plin and Pexp were very similar, but the difference in amplitudes was noticeable (Fig. 4.6 and Fig. 4.7). In eight out of nine subjects, both Plin and Pexp showed a delay during the initial rise to the systolic peak when compared to the corresponding Ptono (Fig. 4.6). This can be explained by the phase lag caused by the viscoelasticity of the 89 4.5. Discussion Ptono Plin Pexp PP (mmHg) 32.6 ± 7.9 30.6± 6.0 33.9 ± 7.4 r = 0.90 r = 0.90 AI −0.28± 0.09 −0.16 ± 0.06 −0.21 ± 0.07 r = 0.71 r = 0.69 TAC (cm3/mmHg) 1.29 ± 0.48 1.18 ± 0.55 1.21 ± 0.56 r = 0.92 r = 0.91 ∆tono−lin ∆tono−exp PP (mmHg) 1.9± 3.6 −1.3± 3.5 p =NS p =NS AI −0.11± 0.06 -0.07 ± 0.06 p < 0.0006 p =NS TAC (cm3/mmHg) 0.11 ± 0.21 −0.01 ± 0.12 p =NS p =NS Table 4.3: Comparisons of pulse pressures, augmentation indices and total arterial compliances obtained in 9 healthy children using applanation tonom- etry and B-mode imaging: Ptono is the carotid artery pressure waveform obtained using applanation tonometry; Plin is the aortic pressure waveform obtained using B-mode imaging with linear calibration; Pexp is the aortic pressure waveform obtained using B-mode imaging with exponential cali- bration; PP is the pulse pressure; AI is the augmentation index; TAC is the total arterial compliance; ∆ indicates the difference; r is the correlation between Ptono and Plin, or Ptono and Pexp. 90 4.5. Discussion 0 0.2 0.4 0.6 0.8 60 65 70 75 80 85 90 95 100 105 P re ss u re ( m m H g ) Time (s) P lin P exp P tono Figure 4.6: The pressure waveforms: Ptono is the carotid artery pressure waveform obtained using applanation tonometry; Plin is the aortic pressure waveform obtained using B-mode imaging with linear calibration; Pexp is the aortic pressure waveform obtained using B-mode imaging with exponen- tial calibration. arterial wall. Both Plin and Pexp also appeared to be “flatter” than the corresponding Ptono. This is particularly obvious in some subjects (Fig. 4.7). Similar observation was obtained in a study comparing the carotid artery distension waveform obtained using echo-tracking and the carotid artery pressure waveform obtained using applanation tonometry [86]. This can again be explained by the phase lag caused by the viscoelasticity and also the nonlinear elasticity. Because of the nonlinear elasticity, the diameter increases by a lesser amount at higher pressure and causes the distension waveform to be “flatter”. The frequency components of Plin and Pexp were only moderately cor- related with that of Ptono. This can be explained by the delay in peak and the “flatter” shape of the distension waveform. AI estimated from Plin and Pexp were also only moderately correlated with that estimated from Ptono. This is likely due to the “flatter” shape of the distension waveform. However, the difference between the AI estimated from Pexp and that esti- mated from Ptono was only 0.07 ± 0.06 and was statistically insignificant. It should also be noted that the inflection point is more difficult to identify for the younger subjects. This is because the reflected wave arrives at the 91 4.5. Discussion 0 0.2 0.4 0.6 0.8 70 75 80 85 90 95 100 105 P re ss u re ( m m H g ) Time (s) P lin P exp P tono Figure 4.7: The pressure waveforms: Ptono is the carotid artery pressure waveform obtained using applanation tonometry; Plin is the aortic pressure waveform obtained using B-mode imaging with linear calibration; Pexp is the aortic pressure waveform obtained using B-mode imaging with exponen- tial calibration; the waveforms estimated from the distension waveforms are “flatter” than that obtained using applanation tonometry. 92 4.5. Discussion ascending aorta during late systole or early diastole, and is “buried” by the incident wave. TAC estimated from Plin and Pexp were highly correlated with that estimated from Ptono. This is because the calculation of TAC involves only the diastolic decay of the pressure waveform, and is not affected by the delay in peak and the “flatter” shape of the distension waveform. The parameters of the semi-automatic edge detection algorithm includ- ing β, θ and ϕ are chosen according to the usual anatomical size of the aortic arch of the children of the particular age group in this study. For adults, a different set of parameters may be needed. The parameter γ of the adaptive block matching scheme is chosen according to the overall image quality of the B-mode images collected in this study. A future study of a large number of subjects and images is necessary to validate the choice of γ. According to the Association for the Advancement of Medical Instru- mentation (AAMI), the criterion for device performance evaluated against the commonly used cuff/stethoscope auscultation method or intra-arterial pressure allows a difference of 5±8 mmHg [174]. PP estimated from Plin and Pexp are well within this acceptability limit if PP estimated from Ptono is used as the reference. Therefore, PP estimated from Plin and Pexp are po- tential surrogates to PP estimated from Ptono. However, it should be noted that the AAMI criterion is for device performance evaluated against either the commonly used cuff/stethoscope auscultation method or intra-arterial pressure, but not applanation tonometry. Kelly and Fitchett [82] proposed a calibration scheme for the carotid artery pressure waveform obtained using applanation tonometry, and the calibrated carotid artery pressure waveform was shown to be an accurate es- timate of the aortic pressure waveform. Since then, this calibration scheme has been popularly used. This calibration scheme is the same as the cal- ibration scheme used for Ptono and Plin in this study, except that MP is calculated as the area integral divided by the duration of the brachial artery pressure waveform. In this study, MP was calculated from the traditional formula (4.6). As the same calibration scheme was used for Ptono and Plin, the error in the measurement of MP (ǫm) does not affect the comparison between Ptono and Plin. However, it is still of considerable interest how ǫm affects the estimation of PP. Since it is assumed that the true aortic pressure waveform is linearly related to the uncalibrated tonometric waveform and the uncalibrated dis- 93 4.6. Conclusion and future work tension waveform, then: SPAO −DP MP−DP = s− d m− d (4.12) where SPAO is the aortic systolic pressure, m is the area integral over the duration of the uncalibrated waveform, s and d are the maximum and minimum of the uncalibrated waveform respectively. SPAO can then be calculated from MP, DP and the uncalibrated waveform using the following formula: SPAO = MD−DP m− d (s− d) + DP (4.13) ǫm propagates towards the estimation of SPAO as follow: SPAO + ǫs = (MP + ǫm)−DP m− d (s− d) + DP (4.14) = SPAO + ǫm m− d (s− d) ǫs = ǫm m− d (s− d) (4.15) where ǫs is the error in the aortic systolic pressure. Assuming that MP is roughly the sum of one third of SP and two thirds of DP, then the term (s-d)/(m-d) is about 3. MP calculated using the traditional formula and that obtained invasively at the brachial artery were compared by Bos et al. [21] and the difference was 5.0 ± 2.3 mmHg. Suppose ǫm is 5 mmHg then ǫs is approximately 15 mmHg. Such error is out of the acceptability limit. However, the calibration scheme proposed by Kelly and Fitchett [82] does not use the invasive MP as the study by Bos et al. [21]. The difference in MP calculated from the traditional formula and that from the area integral of the brachial artery pressure waveform may be smaller than 5 mmHg, if the same sphygmomanometric measurements are used in both calculations. Another non-invasive alternative is to use an automatic oscillometric device that provides direct readout of MP. MP at the brachial artery obtained using such device has been compared to the invasive MP at the aorta [148] and the difference was 0.79 mmHg. Using the oscillometric method to obtain MP may improve the accuracy of the calibration scheme. 4.6 Conclusion and future work A method to estimate the aortic pressure waveform using B-mode imaging is developed. The carotid artery pressure waveform obtained using appla- nation tonometry was used as the reference method. The aortic distension 94 4.6. Conclusion and future work waveform showed a delay during the initial rise to the systolic peak and had a “flatter” shape when compared to the carotid artery pressure waveform. The exponential calibration scheme may be better than the linear calibration scheme since the differences in AI and TAC measured from Ptono and Pexp were smaller and insignificant (Table 4.3). The exponentially-calibrated aor- tic pressure waveform was similar to the tonometric carotid artery pressure waveform. PP and the TAC estimated from the exponentially-calibrated aortic pressure waveform could be used as the surrogate to that estimated from the tonometric carotid artery pressure waveform. AI estimated from the exponentially-calibrated aortic pressure waveform showed moderate ac- curacy. The size of this study was small, and a large-scale study is needed to confirm the validity of the proposed method. The repeatability and repro- ducibility of the proposed method has to be confirmed. The subjects in this study were young, while the augmentation index is often measured in older and hypertensive subjects. It is necessary to investigate the possibility of measuring the augmentation index using the proposed method in older and hypertensive subjects. 95 Chapter 5 Conclusions and Future Work 5.1 Contributions 5.1.1 Model-based estimation of the aortic pulse wave velocity In Chapter 2, a technique was developed to estimate the aortic pulse wave velocity using non-invasive data. The technique is based on a transmission line model and only requires measurements at one site. The model parame- ters are identified using the white-box identification and the prediction-error method. There have been many studies of the arterial system, of human and other mammals, based on the transmission line model with variations in the types of lines, terminations and configurations [15, 24, 25, 28, 44, 78, 117]. However, few studies investigated the subjects of parameter identification. Several studies [24, 25, 28] used the prediction-error method to identify the model parameters. Although the technique presented in Chapter 2 is also based on the prediction-error method, there are several fundamental differ- ences when compared to the techniques proposed in the previous studies. In Chapter 2, the model is formulated in the time domain and its rep- resentation is based on the finite element method. Three parameters, unit line capacitance, load capacitance and line length, are identified using the prediction-error method. As the error function may not be convex, initial values of the parameters are calculated based on the pressure and flow data, and physiological information of the subjects. These initial values are spe- cific to the subject and therefore should be in closer proximity of the true solution, as compared to the initial values chosen from some fixed ranges. The ratio of the line compliance to the load compliance and the pulse wave velocity are calculated from the identified parameters. This technique was tested using non-invasive data obtained from healthy children and children with Marfan syndrome. It was found that the ratio of the line compliance to 96 5.1. Contributions the load compliance was lower in the children with Marfan syndrome than healthy children; and the pulse wave velocity was higher in the children with Marfan syndrome than healthy children. Both observations matched the theoretical expectation. The uniform lossless transmission model used in [25] is the same as that used in Chapter 2, but the pressure and flow waveforms are measured at the descending aorta instead. In [25], the model is formulated in the fre- quency domain and represented by its input impedance. Three parame- ters, the pulse transit time, the characteristic impedance of the line and the time constant of the load impedance are identified using the prediction-error method. It is suggested that the identification problem has infinite solution, and the smallest solution should be chosen as it is the most physiologically relevant. This technique was only tested using invasive data obtained from anesthetized dog, and has yet been applied on human. In [24, 28], a T-tube model is used to represent the arterial system. To identify the parameters, the descending aortic flow has to be estimated. Because of the different identifiable model parameters, the technique presented in Chapter 2 is more clinically useful. Traditionally, pulse wave velocity is estimated using measurements ob- tained at two sites. As presented in Chapter 2, the two-site technique re- quires high accuracy of the measurement of the pulse transit time. In ad- dition, distance measurement can be difficult in some subjects. Recently, a new technique known as Arteriograph has been proposed to measure aor- tic pulse wave velocity [11]. It only requires pressure measurement at the brachial artery. However, the distance measurement from the jugulum to the symphysis is still required. In addition, the theoretical basis of this technique has been challenged [159]. The aortic pulse wave velocity estimation technique presented in Chapter 2 is a one-site technique, and has been shown to be able to differentiate a group of children with Marfan syndrome from a group of healthy children. It is based on the transmission line model, which has been regarded as a suitable model for the arterial system [94]. It is also a non-invasive technique, making it potentially useful in routine clinical practice. 5.1.2 Minimally-invasive monitoring of the cardiac output In chapter 3, a technique for continuous cardiac output monitoring using radial artery pressure was developed. The technique uses a state-space representation of a transmission line model to represent the arterial sys- tem. The model parameters are identified using a novel grey-box identifica- 97 5.1. Contributions tion method. The identifiability of the proposed state-space representation was investigated. The performance of the proposed grey-box identification method was tested using a Gaussian input and a simulated aortic waveform input. The proposed state-space representation is based on the finite element method [96]. Compared to the transfer function of the transmission line model, the state-space representation has more flexibility. In the state-space representation, the model parameters, including total line compliance, total line inductance, load compliance and load resistance, are separate entries of the system matrix. In contrast, the transfer function is a lumped represen- tation of these parameters. This feature of the state-space representation is especially useful during monitoring of the cardiac output. For example, if additional information of the compliance is available during monitoring, the model can be adjusted accordingly because the initial values of the param- eters are known. In the proposed grey-box identification, subspace identification is ap- plied to identify an empirical model, then the physical model parameters are identified from the empirical model through iterative bilinear regres- sions. The proposed method is based on the framework by Xie and Ljung [175]. In Chapter 3, it was showed how this framework can be applied when the system matrices have unknown entries that are linear functions of a small number of physical parameters. The proposed method is more computational efficient than the white-box method presented in Chapter 2. In the presented white-box method, an optimization algorithm, such as the gradient-descent method and the Newton conjugate-gradient method, has to be used to identify the model parameters. In most cases, these optimiza- tion algorithms are more computationally expensive than a series of bilinear regressions. The proposed technique was applied to 5 infants during their post- operative stay in the Intensive care unit. The accuracy of this technique was shown to be clinically acceptable during a 5 hours re-calibration-free period. As presented in Section1.3.3, PulseCO and PiCCO are the current meth- ods that monitor the cardiac output using the radial artery pressure wave- form. Both of these techniques require invasive calibration using thermodilu- tion or transpulmonary thermodilution. In contrast, for the proposed tech- nique, calibration is performed non-invasively using transthoracic Doppler ultrasound. Recently, techniques based on multichannel blind system identification have been proposed to monitor cardiac output using multiple peripheral 98 5.1. Contributions pressure waveforms [107,155]. Compared to PulseCO and PiCCO, these techniques are based on more sophisticated identification method. However, calibration is still required by these techniques, and an additional peripheral pressure measurement is also needed. Up until now, these techniques are yet to be proved to be superior to PulseCO and PiCCO for cardiac output monitoring in human. 5.1.3 Estimation of the aortic pressure waveform using B-mode imaging In Chapter 4, a method to measure aortic pressure waveform using B-mode imaging was presented. Recently, the idea of estimating the carotid artery pressure waveform from the carotid artery distension waveform has been explored by several studies [86,161,166]. In these studies, echo-tracking was used to measure the carotid artery distension waveform. Given the high spatial resolution that can be achieved by echo-tracking, the carotid artery distension waveform can be measured reliably. In contrast, the reproducibil- ity and repeatability of the edge detection algorithm presented in chapter 4 has yet to be investigated. The method presented in chapter 4 may have an advantage over the echo-tracking method as it measures the aortic pressure waveform using the distension waveform obtained at the aorta. The carotid artery pressure waveform may not be a suitable surrogate as the aortic pressure waveform in the case of carotid stenosis. The echo-tracking method has yet to be used to measure the aortic pressure waveform. 5.1.4 Suggestions for future work It has been shown that the model-based pulse wave velocity estimation method proposed in this thesis can differentiate a group of children with Marfan syndrome from a group of normal children. However, since the true aortic stiffness is not known, the accuracy of this method is yet to be con- firmed. The true aortic stiffness can be obtained by measuring the change in the aortic diameter and the corresponding change in the pressure at the same location invasively. A follow-up animal study in which the true aortic stiffness is altered and the correlation between the true aortic stiffness and the model-based pulse wave velocity is found may confirm the accuracy of the proposed method. For the cardiac output monitoring study, a major limitation is that the number of patients participated is small. Studies with a much larger sam- 99 5.1. Contributions ple size are required to further validate this technique. In addition, it is useful to determine for how long this technique can maintain an acceptable accuracy without re-calibration. Further studies to acquire multiple CO determinations at different time points in a longer time intervals from each patient would be valuable. This technique should also be validated on differ- ent patient populations and under different conditions (e.g. hemodynamic instability). Using the proposed grey-box identification method, the patient-specific transmission line model is identified and the aortic flow can be estimated from the radial artery pressure. Aortic pressure may also be estimated from the identified model and the radial artery pressure. A non-invasive method to measure the aortic pressure waveform is desirable, and the generalized transfer function method has been used to estimate the aortic pressure from the radial pressure. The possibility of using the patient-specific transmission line model instead of the generalized transfer function should be investigated. The method for aortic pressure waveform measurement presented in this thesis is strongly dependent on the accuracy of the method for the aortic distension waveform measurement. The reproducibility, the repeatability, the inter- and intra- observer errors of the method for the aortic distension waveform measurement should be further investigated. 100 References 1. E. Agabiti-Rosei, G. Mancia, M. F. O’Rourke, M. J. Roman, M. E. Safar, H. Smulyan, J.-G. Wang, I. B. Wilkinson, B. Williams and C. Vlachopoulos, “Central blood pressure measurements and antihy- pertensive therapy: a consensus document,” Hypertension, vol. 50, pp.154-160, 2007. 2. E. L. Alderman, A. Branzi, W. Sanders and B. W. Brown, “Evaluation of the pulse-contour method of determining stroke volume in man,” Circulation, vol. 46, pp. 546-558, 1972. 3. R. Armentano, J. L. Megnien, A. Simon, F. Bellenfant, J. Barra and J. Levenson, “Effects of hypertension on viscoelasticity of carotid and femoral arteries in humans,” Hypertension, vol. 26, pp. 48-54, 1995. 4. R. G. Asmar, P. C. Brunel, B. M. Pannier, P. J. Lacolley and M. E. Safar, “Arterial distensibility and ambulatory blood pressure monitor- ing in essential hypertension,” Am J Cardiol, vol. 61, pp. 1066-1070, 1988. 5. R. Asmar, A. Benetos, G. London, C. Hugue, Y. Weiss, J. Topouch- ian, B. Laloux and M. Safar, “Aortic distensibility in normotensive, untreated and treated hypertensive patients,” Blood Press, vol. 4, pp. 48-54, 1995. 6. R. Asmar, A. Benetos, J. Topouchian, P. Laurent, B. Pannier, A.-M. Brisac, R. Target and B. I. Levy, “Assessment of arterial distensibility by automatic pulse wave velocity measurement: validation and clinical application studies,” Hypertension, vol. 26, pp. 485-490, 1995. 7. R. G. Asmar, G. M. London, M. E. O’Rourke and M. E. Safar, “Im- provement in blood pressure, arterial stiffness and wave reflections with a very-low-dose Perindopril/Indapamide combination in hyper- tensive patient,” Hypertension, vol. 38, pp. 922-926, 2001. 101 References 8. A. P. Avolio, S.-G. Chen, R.-P. Wang, C.-L. Zhang, M.-F. Li and M. F. O’Rourke, “Effects of aging on changing arterial compliance and left ventricular load in a northern Chinese urban community,” Circulation, vol. 68, pp. 50-58, 1983. 9. A. P. Avolio, F.-Q. Deng, W.-Q. Li, Y.-F. Luo, Z.-D. Huang, L.-F. Xing and M. F. O’Rourke, “Effects of aging on arterial distensibility in populations with high and low prevalence of hypertension: comparison between urban and rural communities in China,” Circulation, vol. 71, pp. 202-210, 1985. 10. J.-P. Baguet, S. Douchin, H. Pierre, A.-M. Rossignol, M. Bost and J.- M. Mallion, “Structural and functional abnormalities of large arteries in the Turner syndrome,” Heart, vol. 91, pp. 1442-1446, 2005. 11. J. Baulmann, U. Schillings, S. Rickert, S. Uen, R. Düsing, M. Illyes, A. Cziraki, G. Nickenig and T. Mengden, “A new oscillometric method for assessment of arterial stiffness: comparison with tonometric and piezo-electronic methods,” J Hypertens, vol. 26, pp. 523-528, 2008. 12. D. H. Bergel, “The static elastic properties of the arterial wall,” J Physiol, vol. 156, pp. 445-457, 1961. 13. D. H. Bergel, “The dynamic elastic properties of the arterial wall,” J Physiol, vol. 156, pp. 458-469, 1961. 14. D. S. Berger, J. K. Li, W. K. Laskey and A. Noordergraaf, “Repeated reflection of waves in the systemic arterial system,” Am J Physiol, vol. 264, pp. 269-291, 1993. 15. D. S. Berger, J. K. Li and A. Noordergraaf, “Differential effects of wave reflections and peripheral resistance on aortic blood pressure: a model-based study,” Am J Physiol, vol. 266, pp. 1626-1642, 1994. 16. G. Biancofiore, L. A. H. Critchley, A. Lee, L. Bindi, M. Bisà, M. Es- posito, L. Meacci, R. Mozzo, P. DeSimone, L. Urbani and F. Filipponi, “Evaluation of an uncalibrated arterial pulse contour cardiac output monitoring system in cirrhotic patients undergoing liver surgery,” Br J Anaesth, vol. 102, pp. 47-54, 2009. 17. J. Blacher, R. Asmar, S. Djane, G. M. London and M. E. Safar, “Aortic pulse wave velocity as a marker of cardiovascular risk in hypertensive patients,” Hypertension, vol. 33, pp. 1111-1117, 1999. 102 References 18. R. D. Bland, W. C. Shoemaker, E. Abraham, J. C. Cobo, “Hemody- namic and oxygen transport patterns in surviving and nonsurviving postoperative patients,” Crit Care Med, vol.13, pp. 85-90, 1985. 19. P. Boutouyrie, A. I. Tropeano, R. Asmar, I. Gautier, A. Benetos, P. Lacolley, S. Laurent, “Aortic stiffness is an independent predictor of primary coronary events in hypertensive patients: a longitudinal study,” Hypertension, vol. 39, pp. 10-15, 2002. 20. M. Boyle, L. Steel, G. M. Flynn, M. Murgo, L. Nicholson, M. O’Brien and D. Bihari, “Assessment of the clinical utility of an ultrasonic mon- itor of cardiac output (the USCOM) and agreement with thermodi- lution measurement,” Crit Care Resusc, vol. 11, pp. 198-203, 2009. 21. W. Bos, E. Verrij, H. Vincent, B. Westerhof, G. Parati and G. van Montfrans G, “How to assess mean blood pressure properly at the brachial artery level,” J Hypertens, vol. 25, pp. 751-755, 2007. 22. T. J. Bradley, J. E. Potts, M. T. Potts, A. M. DeSouza and G. G. S. Sandor, “Echocardiographic Doppler assessment of the biophysi- cal properties of the aorta in pediatric patients with the Marfan syn- drome,” Am J Cardiol, vol. 96, pp. 1317-1321, 2005. 23. J. C. Bramwell and A. V. Hill, “Velocity of transmission of the pulse- wave,” Lancet, vol. 199, pp. 891-892, 1922. 24. Burattini, R., and K. B. Campbell. Modified asymmetric T-tube model to infer arterial wave reflection at the aortic root. IEEE Trans Biomed Engg, vofl. 36, pp. 805-814, 1989. 25. R. Burattini K. B. Campbell, “Physiological relevance of uniform elas- tic tube-models to infer descending aortic wave reflection: a problem of identifiability,” Ann Biomed Eng, vol. 28, pp. 512-523, 2000. 26. R. Burattini and P. O. Di Salvia, “Development of systemic arterial mechanical properties from infancy to adulthood interpreted by four- element windkessel models,” J Appl Physiol, vol. 103, pp.66–79, 2007. 27. K. B. Campbell, L. C. Lee, H. F. Frasch and A. Noordergraaf, “Pulse reflection sites and effective length of the arterial system,” Am J Phyiol, vol. 256, pp. 1684-1689, 1989. 103 References 28. K. B. Campbell, R. Burattini, D. L. Bell, R. D. Kirkpatrick and G. G. Knowlen, “Time-domain formulation of asymmetric T-tube model of arterial system,” Am J Physiol, vol. 258, pp. 1761-1774, 1990. 29. M. Cecconi, J. Fawcett, R. M. Grounds and A. Rhodes, “A prospective study to evaluate the accuracy of pulse power analysis to monitor cardiac output in critically ill patients,” BMC Anesth, vol. 8, 2008. 30. J. S Chan, D. Segara and Priya Nair, “Measurement of cardiac output with a non-invasive continuous wave Doppler device versus the pul- monary artery catheter: a comparative study,” Crit Care Resusc, vol. 8, pp, 309-314, 2006. 31. R. Chand, Y. Mehta and N. Trehan, “Cardiac output estimation with a new Doppler device after off-pump coronary artery bypass surgery,” J Cardiothorac Vasc Anesth, vol. 20, pp. 315-319, 2006. 32. C.-H. Chen, E. Nevo, B. Fetics, P. H. Pak, F. C. P. Yin, L. Maughan and D. A. Kass, “Estimation of central aortic pressure waveform by mathematical transformation of radial tonometry pressure,” Circula- tion, vol. 95, pp. 1827-1836, 1997. 33. G. C. Cloud, C.Rajkumar, J. Kooner, J. Cooke and C. J. Bulpitt, “Estimation of central aortic pressure by SphygmoCor requires intra- arterial peripheral pressures,” Clin Sci, vol. 105, pp. 219-225, 2003. 34. A. F. Connors Jr, T. Speroff, N. V. Dawson, C. Thomas, F. E. Harrell Jr, D. Wagner, N. Desbiens, L. Goldman, A. W. Wu, R. M. Califf, W. J. Fulkerson Jr, H. Vidaillet, S. Broste, P. Bellamy, J. Lynn and W. A. Knaus, “The effectiveness of right heart catheterization in the initial care of critically ill patients,” JAMA, vol. 276, pp. 889-897, 1996. 35. A. Corley, A. G. Barnett, D. Mullany and J. F. Fraser, “Nurse-determined assessment of cardiac output. Comparing a non-invasive cardiac out- put device and pulmonary artery catheter: a prospective observational study,” Int J Nurs Stud, vol. 46, pp. 1291-1297, 2009. 36. L. A. H. Critchley and J. A. J. H. Critchley, “A meta-analysis of studies using bias and precision statistics to compare cardiac output estimation techniques,” J Clin Monit, vol. 15, pp.85–91, 1999. 104 References 37. J. I. Davies, M. M. Band, S. Pringle, S. Ogston and A. D. Struthers, “Peripheral blood pressure measurement is as good as applanation tonometry at predicting ascending aortic blood pressure,” J Hypertens, vol. 21, pp. 571-576, 2003. 38. J. De Backer, B. Loeys and A. De Paepe, “Marfan and Marfan-like syndromes,” Artery Res, vol. 3, pp. 9-16, 2009. 39. K. De Cock and B. De Moor, “Subspace identification methods,” Con- trol systems, robotics and automation, EOLSS, UNESCO encyclopedia of life support systems, vol. 1, pp. 933–979, 2003. 40. R. H. Dietzman, R. A. Ersek, C. W. Lillehei, A. R. Castaneda and R. C. Lillehei, “Low output syndrome: recognition and treatment,” J Thorac Cardiovasc Surg, vol. 57, pp. 138-150, 1969. 41. J. A. Dineley, N. McDicken and P. R. Hoskins, “Factors affectng the ar- terial distension waveform derived from tissue doppler imaging (TDI): an in vintro study on precision,” Ultrasound Med Biol, vol.33, pp. 1123-1131, 2007. 42. H. G. M. Dötsch and P. M. J. Van den Hof, “Test for local struc- ture identifiability of high-order non-linearly parameterized state space models,” Automatica, vol. 32, pp. 873–883, 1996. 43. G. M. Drzewiecki, J. Melbin and A. Noordergraaf, “Arterial tonom- etry: review amnd analysis,” J Biomechanics, vol. 16, pp. 141-152, 1983. 44. Einav, S., S. Aharoni, and M. Manoach. Exponentially tapered trans- mission line model of the arterial system. IEEE Trans. Biomed. Eng. 35:333-339, 1988. 45. D. J. Eveson, T. G. Robinson, N. S. Shan, R. B. Panerai, S. K. Paul and J. F. Potter, “Abnormalities in cardiac baroreceptor sensitivity in acute ischaemic stroke patients are related to aortic stiffness,” Clin Sci, vol. 108, pp. 441-447, 2005. 46. B. Fetics, E. Nevo, C.-H. Chen and D. A. Kass, “Parameteric model derivation of transfer function for noninvasive estimation of aortic pres- sure by radial tonometry,” IEEE Trans Biomed Eng, vol. 46, pp. 698- 706, 1999. 105 References 47. U. Fakler, Ch. Pauli, G. Balling, H. P. Lorenz, A. Eicken, M. Hennig, and J. Hess, “Cardiac index monitoring by pulse contour analysis and thermodilution after pediatric cardiac surgery,” J Thorac Cardiovasc Surg, vol.133, pp. 224–228, 2007. 48. T. W. Felbinger, D. A. Reuter, H. K. Eltzschig, K. Moerstedt, O. Goedje and A. E. Goetz, “Comparison of pulmonary arterial ther- modilution and arterial pulse Contour analysis: evaluation of a new algorithm,” J Clin Anesth, vol. 14, pp. 296-301, 2002. 49. T. W. Felbinger, D. A. Reuter, H. K. Eltzschig, J. Bayerlein and A. E. Goetz, “Cardiac index measurements during rapid preload changes: a comparison of pulmonary artery thermodilution with arterial pulse contour analysis,” J Clin Anesth, vol. 17, pp. 241-248, 2005. 50. E. A. Fisher, J. A. Stahl, J. H. Budd, M. E. Goldman, “Transesophageal echocardiography: procedures and clinical application,” JACC, vol. 18, pp. 1333-1348, 1991. 51. J. A. Franciosa, “Evaluation of the CO2 rebreathing cardiac out- put method in seriously ill patients,” Circulation, vol. 55, pp. 449- 455,1977. 52. D. N. Ghista, G. Jayaraman and H. Sandler, “Analysis for the non- invasive determination of arterial properties and for the transcuta- neous continuous monitoring of arterial blood pressure,”Med Biol Eng Comput, vol. 16, pp. 715-726, 1979. 53. O. Gödje, C. Thiel, P. Lamm, H. Reichenspurner, C. Schmitz, A. Schütz and B. Reichart, “Less invasive, continuous hemodynamic mon- itoring during minimally invasive coronary surgery,” Ann Thorac Surg, vol. 68, pp. 1532-1536, 1999. 54. O. Gödje, K. Höke, A. E. Goetz, T. W. Felbinger, D. A. Reuter, B. Reichart, R. Friedl, A. Hannekum, and U. J. Pfeiffer, “Reliability of a new algorithm for continuous cardiac output determination by pulse- contour analysis during hemodynamic instability,” Crit Care Med, vol. 30, pp. 52–59, 2002. 55. S. Golemati, A. Sassano, M. J. Lever, A. A. Bharath, S. Dhanjiu and A. N. Nicolaides, “Carotid artery wall motion estimated from B-mode ultrasound using region tracking and block matching,” Ultrasound Med Biol, vol. 29, pp. 387-399, 2003. 106 References 56. S. Graf, J. Gariepy, M. Massonneau, R. L. Armentano, S. Mansour, J. G. Barra, A. Simon and J. Levenson, “Experimental and clinical validation of arterial diameter waveform and intmal media thickness obtained from B-mode ultrasound image processing,” Ultrasound Med Biol, vol. 25, pp. 1353-1363, 1999. 57. J. C. Greenfield and D. J. Patel, “Relation between pressure and di- ameter in the ascending aorta of man,” Circ Res, vol. 10, pp. 778-781, 1962. 58. M. Groenink, A. de Roos, B. J. M. Mulder, J. A. E. Spaan an E. E. van der Wall, “Changes in aortic distensibility and pulse wave velocity assessed with magnetic resonance imaging following beta-blocker ther- apy in the Marfan syndrome,” Am J Cardiol, vol. 82, pp. 203-208, 1998. 59. M. Groenink, A. De Roos, B. J. M. Mulder, B. Verbeeten Jr., J. Tim- mermans, A. H. Zwinderman, J. A. E. Spaan and E. E. van der Wall, ” Biophysical properties of the normal-sized aorta in patients with Mar- fan Syndrome: evaluation with MR flow mapping,” Radiology, vol. 219, pp. 535-540, 2001. 60. J.-O. Hahn, A. T. Reisner, and H. H. Asada, “Blind identification of two-channel IIR systems with application to central cardiovascular monitoring,” J-DSMC, vol. 131, 2009. 61. B. A. Haluska, L. Jeffriess, P. M. Mottram, S. G. Carlier and T. H. Marwick, “A new technique for assessing arterial pressure wave forms and central pressure with tissue Doppler,” Cardiovasc Ultrasound, vol. 5, pp. 6, 2007. 62. W. F. Hamilton and J. W. Remington, “The measurement of the stroke volume from the pressure pulse,” Am J Physiol, vol. 148, pp. 14-24, 1947. 63. T. T. Hamilton, L. M. Huber, and M. E. Jessen, “PulseCO: A less- invasive method to monitor cardiac output from arterial pressure after cardiac surgery,” Ann Thorac Surg, vol. 74, pp. 1408–1412, 2002. 64. O. Hamzaoui, X. Monnet, C. Richard, D. Osman, D. Chemla, and J.-L. Teboul, “Effects of changes in vascular tone on the agreement between 107 References pulse contour and transpulmonary thermodilution cardiac output mea- surements within an up to 6-hour calibration-free period,” Crit Care Med, vol. 36, pp. 434–440, 2008. 65. C. S. Hayward, D. C. Knight, B. G. Wren and R. P. Kelly, “Effect of hormone replacement therapy on non-invasive cardiovascular haemo- dynamics,” J Hypertens, vol. 15, pp. 987-993, 1997. 66. A. P. G. Hoeks, P. J. Brands, F. A. M. Smeets and R. S. Reneman, “Assessment of the distensibility of superficial arteries,” Ultrasound Med Biol vol. 16, pp. 121-128, 1990. 67. C. K. Hofer, M. Cecconi, G. Marx, and G. della Rocca, “Minimally invasive haemodynamic monitoring,” Eur. J. Anaesthesiol, vol. 26, pp. 996–1002, 2009. 68. S.A. Hope, D. B. Tay, I. T. Meredith and J. D. Cameron, “Use of arterial transfer functions for the derivation of aortic waveform char- acteristics,” J Hypertesn, vol. 21, pp. 1299-1305, 2003. 69. S. A. Hope, I. T. Meredith, D. Tay and J. D. Cameron, “ ‘Generaliz- ability’ of a radial-aortic transfer function for the derivation of central aortic waveform parameters,” J Hypertens, vol. 25, pp. 1812-1820, 2007. 70. S. A. Hope, I. T. Meredith and J. D. Cameron, “Arterial transfer functions and the reconstruction of central aortic waveforms: myths, controversies and misconceptions,” J Hypertens, vol. 26, pp. 4-7, 2008. 71. I. G. Horváth, A. Németh, Z. Lenkey, N. Alessandri, F. Tufano, P. Kis, B. Gaszner and A. Cziráki, ”Invasive validation of a new oscillomet- ric device (Arteriograph) for measuring augmentation index, central blood pressure and aortic pulse wave velocity,” J Hypertens, vol. 28, 2010. 72. T. Imura, K. Yamamoto, T. Satoh, K. Kanamori, T. Mikami and H. Yasuda, “In vivo viscoelastic behavior in the human aorta,” Circ Res, vol. 66, pp. 1413-1419, 1990. 73. R. N. Isnard, B. M. Pannier, S. Laurent, G. M. London, B. Diebold and M. E. Safar, “Pulsatile diameter and elastic modulus of the aortic arch in essential hypertension: a noninvasive study,” JACC, vol. 13, pp. 399-405, 1989. 108 References 74. J. R. C. Jansen, K. H. Wesseling, J. J. Settels and J. J. Schreuder, “Continuous cardiac output monitoring by pulse contour during car- diac surgery,” Eur Heart J, vol. 11(Suppl. 1), pp. 26-32, 1990. 75. N. A. Jatoi, A. Mahmud, K. Bennett and J. Feely, “Assessment of arterial stiffness in hypertension: comparison of oscillometric (Arte- riograph), piezoelectronic (Complior) and tonometric (SphygmoCor) techniques,” J Hypertens, vol. 27, pp. 2186-2191, 2009. 76. G. Jondeau, P. Boutouyrie, P. Lacolley, B. Laloux, O. Dubourg, J. Bourdarias and S. Laurent, “Central pulse pressure is a major deter- minant of ascending aorta dilation in Marfan syndrome,” Circulation, vol. 99, pp. 2677-2681, 1999. 77. M. Karamanoglu, M. F. O’Rourke, A. P. Avolio and R. P. Kelly, “An analysis of the relationship between aortic and peripheral upper limb pressure waves in man,” Eur Heart J, vol. 14, pp. 160-167, 1993. 78. M. Karamanoglu, D. E. Gallagher, A. P. Avolio and M. F. O’Rourke, “Functional origin of reflected pressure waves in a multibranched model of the human arterial system,” Am J Physiol, vol. 267, pp. 1681-1688, 1994. 79. M. Karamanoglu, D. E. Gallagher, A. P. Avolio and M. F. O’Rourke, “Pressure wave propagation in a multibranched model of the human of the human upper limb,” Am J Physiol, vol. 269, pp. 1363-1369, 1995. 80. R. Kelly, C. Hayward, J. Ganis, J. Daley, A. Avolio and M. O’Rourke, “Noninvasive registration of the arterial pressure pulse waveform using high-fidelity applanation tonometry,” J Vasc Med Biol, vol. 1, pp. 142-149, 1989. 81. R. Kelly, M. Karamanoglu, A. Avolio and M. O’Rourke, “Noninvasive carotid pressure wave registration as an indicator of ascending aortic pressure,” J Vasc Med Biol, vol. 1, pp. 241-241, 1989. 82. R. Kelly and D. Fitchett, ”Noninvasive determination of aortic input impedance and external left ventricular power output: a validation and repeatability study of a new technique,” J Am Coll Cardiol, vol. 20, pp. 952-963, 1992. 109 References 83. R. P. Kelly, S. C. Millasseau, J. M. Ritter and P. J. Chowienczyk, “Vasoactive drugs influence aortic augmentation index independently of pulse-wave velocity in healthy men,” Hypertension, vol. 37, pp. 1429-1433, 2001. 84. J. J. Kim, W. J. Dreyer, A. C. Chang, J. P. Breinholt III, and R. G. Grifka, “Arterial pulse wave analysis: An accurate means of determin- ing cardiac output in children,” Pediatr Crit Care Med, vol. 7, pp. 532–535, 2006. 85. T. Kind, T. J. C. Faes, J.-W. Lankhaar, A. Vonk-Noordegraaf, and M. Verhaegen, “Estimation of three- and four- element Windkessel parameters using subspace model identification,” IEEE Trans Biomed Eng, vol. 57, pp. 1531–1538, 2010 86. J. Kips, F. Vanmolkot, D. Mahieu, S. Vermeersch, I. Fabry, J. de Hoon, L. Van Bortel and P. Segers, “The use of diameter distension waveforms as an alternative for tonometric pressure to assess carotid blood pressure,” Physiol Meas, vol. 31, pp 543-553, 2010. 87. W. Knirsch, O. Kretschmar, M. Tomaske, K. Stutz, N. Nagdyman, C. Balmer, A. Schmitz, D. Béttex, F. Berger, U. Bauersfeld and M. Weiss, “Cardiac output measurement in children: comparison of the ultrasound cardiac output monitor with thermodilution cardiac output measurement,” Intensive Care Med, vol. 34, pp. 1060-1064, 2008. 88. M. J. F. Kool, T. van Merode, R. s. Reneman, A. P. G. Hoeks, H. A. J. Struyker Boudier and L. M. A. B. Van Bortel, “Evaluation of repro- ducibility of a vessel wall movement detector system for assessment of large artery properties,” Cardiovasc Res, vol. 28, pp. 610-614, 1994. 89. E. J. Kroeker and E. H. Wood, “Comparison of simultaneously recorded central and peripheral arterial pressure pulses during rest, exercise and tilted position in man,” Circ Res, vol. 3, pp. 623-632, 1955. 90. P. Lacolley, S. Gautier, O. Poirier, B. Pannier, F. Cambien and A. Benetos, “Nitric oxide synthase gene polymorphisms, blood pressure and aortic stiffness in normotensive and hypertensive subjects,” J Hy- pertens, vol. 16, pp. 31-35, 1998. 91. G. J. Langewouters, K. H. Wesseling and W. J. Goedhard, “The static elastic properties of 45 human thoracic and 20 abdominal aortas in 110 References vitro and the parameters of a new model,” J Biomech, vol. 17, pp. 423-435, 1984. 92. S. Laurent, P. Boutouyrie, R. Asmar,I. Gautier, B. Laloux, L. Guize, P. Ducimetiere and A. Benetos, “Aortic stiffness is an independent predicotr of all-cause and cardiovasular mortality in hypertensive pa- tients,” Hypertension, vol. 37, pp. 1236-1241, 2001. 93. S. Laurent, S. Katsahian, C. Fassot, A.-I. Tropeano I. Gautier, B. Laloux and P. Boutouyrie, “Aortic stiffness is an independent predictor of fatal stroke in essential hypertension,” Stroke, vol. 34, pp. 1203- 1206, 2003. 94. S. Laurent, J. Cockcroft, L. Van Bortel, P. Boutouyrie, C. Giannat- tasio, D. Hayoz, B. Pannier, C. Vlachopoulos, I. Wilkinson and H. Struijker-Boudier, “Expert consensus document on arterial stiffness: methodological issues and clinical applications,” Eur Heart J, vol. 27, pp. 2588-2605, 2007. 95. B. M. Learoyd and M. G. Taylor, “Alterations with age in the vis- coelastic properties of human arterial walls,” Circ Res, vol. 18, pp. 278-292, 1966. 96. S. Y. Lee, A. Konrad, and R. Saldanha, “Lossy transmission line tran- sient analysis by the finite element method,” IEEE Trans Magn, vol. 29, pp. 1730–1732, 1993. 97. N. W. F. Linton and R. A. F. Linton, “Estimation of changes in cardiac output from the arterial blood pressure waveform in the upper limb,” Br J Anaesth, vol. 86, pp. 486-496, 2001. 98. Z. Liu, K. P. Brin and F. C. P. Yin, “Estimation of total arterial compliance: an improved method and evaluation of current methods,” Am J Physiol, vol. 251, pp. 588-600, 1986. 99. L.Ljung, System Identification: Theory for the User. Upper Saddle River, NJ: Prentice-Hall, 1999. 100. G. M. London S. J. Marchais, M. E. Safar, A. F. Genest, A. P. Guerin, F. Metivier, K. Chedid and A. M. London, “Aortic and large artery compliance in end-stage renal failure,” Kidney Int, vol. 37, pp. 137- 142, 1990. 111 References 101. G. M. London, A. P. Guerin, S. J. Marchais, B. Pannier, M. E. Safar, M. Day and F. Metivier, “Cardiac and arterial interactions in end- stage renal disease,” Kidney Int, vol. 50, pp. 600-608, 1996. 102. C. Lyzell, M. Enqvist and L. Ljung, “Handling certain structure in- formation in subspace identification,” in 15th IFAC Symposium on System Identification, Saint-Malo, France, 2009, pp. 90-95. 103. S. Madhavan, W. L. Ooi, H. Cohen and M. H. Alderman, “Relation of pulse pressure and blood pressure reduction to the incidence of myocardial infarction,” Hypertension, vol. 23, pp. 395-401, 1994. 104. A. Mahmud and J. Feely, “Acute effect of caffeine on arterial stiffness and aortic pressure waveform,” Hypertension, vol. 38, pp. 227-231, 2001. 105. G. R. Manecke and W. R. Auger, “Cardiac output determination from the arterial pressure wave: clinical testing of a novel algorithm that does not require calibration,” J Cardiothorac Vasc Anesth, vol. 21, pp. 3-7, 2007. 106. J. Mayer, J. Boldt, M. W. Wolf, J. Lang and S. Suttner, “Cardiac output derived from arterial pressure waveform analysis in patients undergoing cardiac surgery: validity of a second generation device,” Anesth Analg, vol. 106, pp. 867-872, 2008. 107. D. B. McCombie, A. T. Reisner, and H. H. Asada, “Laguerre-model blind system identification: cardiovascular dynamics estimated from multiple peripheral circulatory signals,” IEEE Trans Biomed Eng, vol. 52, pp. 1889–1901, 2005. 108. J. V. McCoy, S. M. Hollenberg, R. P. Dellinger, R. C. Arnold, L. Ru- oss, V. Lotano, P. Peters, J. E. Parrillo and S. Trzeciak, “Continuous cardiac index monitoring: a prospective observational study of agree- ment between a pulmonary artery catheter and a calibrated minimally invasive technique,” Resuscitation, vol. 80, pp. 893-897, 2009. 109. S. Meaume, A. Benetos, O. F. Henry, A. Rudnichi and M. E. Safar, “Aortic pulse wave velocity predicts cardiovascular mortality in sub- jects >70 years of age,” Arterioscler Thromb Vasc Biol, vol. 21, pp. 2046-2050, 2001. 112 References 110. J. M. Meinders and A. P. G. Hoeks, “Simultaneous assessment of di- ameter and pressure waveforms in the carotid artery,” Ultrasound Med Biol, vol. 30, pp. 147-154, 2004. 111. S. C. Millasseau, R. P. Kelly, J. M. Ritter and P. J. Chowienczyk, “De- termination of age-related increases in large artery stiffness by digital pulse contour analysis,” Clin Sci, vol. 103, pp. 371-377, 2002. 112. S. C. Millasseau, A. D. Stewart, S. J. Patel, S. R. Redwood and P. J. Chowienczyk, “Evaluation of carotid-femoral pulse wave velocity; influence of timing algorithm and heart rate,” Hypertension, vol. 45, pp. 222-226, 2005. 113. G. F. Mitchell, L. A. Moye, E. Braunwald, J.-L. Rouleau, V. Bernstein, E. M. Geltman, G. C. Flaker and M. A. Pfeffer, “Sphygmomanometri- cally determined pulse pressure is a powerful independent predictor of recurrent events after myocardial infarction in patients with impaired left ventricular function,” Circulation, vo. 9, pp. 4254-4260, 1997. 114. G. F. Mitchell, H. Parise, E. J. Benjamin, M. G. Larson, M. J. Keyes, J. A. Vita, R. S. Vasan and D. Levy, “Changes in arterial stiffness and wave reflection with advancing age in healthy men and women: the framingham heart study,” Hypertension, vol. 43, pp. 1239-1245, 2004. 115. R. Mukkamala, A. T. Reisner, H. M. Hojman, R. G. Mark, and R. J. Cohen, “Continuous caridac output monitoring by peripheral blood pressure waveform analysis,” IEEE Trans Biomed Eng, vol. 53, pp. 59–67, 2006. 116. J. P. Murgo, N. Westerhof, J. P. Giolma and S. A. Altobelli, “Aortic input impedance in normal man: relationship to pressure wave forms,” Circulation, vol. 62, pp. 105-116, 1980. 117. Myers, L. J., and W. L. Capper. Exponential taper in arteries: an ex- act solution of its effect on blood flow velocity waveforms and impedance. Med. Eng. Phys. 26:147-155, 2004. 118. J. Nemcsik, J. Egresits, T. E. H. Othmane, B C. Fekete, E. Fodor, T. Szabó, Z. Járai, C. Jekkel, I. Kiss and A. Tislér, “Validation of Arteriograph - a new oscillometric device to measure arterial stiffness in patients on maintenance hemodialysis,” Kidney Blood Press Res, vol. 32, pp. 223-229, 2009. 113 References 119. W. W. Nichols and M. F. O’Rourke, McDonald’s Blood Flow in Arter- ies: Theoretical, Experimental and Clinical Principles, London: Hod- der Arnold, 2005. 120. J. Nocedal and S. J. Wright, Numerical Optimization, New York : Springer, 1999. 121. A. Noordergraaf, Circulatory System Dynamics, New York: Academic Press, 1978. 122. M. F. O’Rourke, “Mechanical principles in arterial disease,” Hyper- tension, vol. 26, pp.2-9, 1995. 123. M. F. O’Rourke, J. A. Staessen, C. Vlachopoulos, D. Duprez and G. E. Plante, “Clinical applications of arterial stiffness: definitions and reference values,” Am J Hypertens, vol. 15, pp. 426-444, 2002. 124. B. Pannier, P. Brunel and W. El Aroussy, “Pulse pressure and echocar- diographic findings in essential hypertension,” J Hypertens, vol. 7, pp. 127-132, 1989. 125. N. Patel, M. Dodsworth and J. F. Mill, “Cardiac output measure- ment in newborn infants using the ultrasonic cardiac output infants using the ultrasonic cardiac output monitor: an assessment of agree- ment with conventional echocardiography, repeatability and new user experience,” Arch Dis Child Fetal Neonatal Ed, 2010. 126. A. L. Pauca, M. F. O’Rourke and N. D. Kon, “Propective evaluation of a method for estimating ascending aortic pressure from the radial artery pressure waveform, ” Hypertension, vol. 38, pp. 932-937, 2001. 127. A. C. Perrino Jr, “Cardiac output monitoring by echocardiography: should we pass on Swan-Ganz catheters,” Yale J Biol Med, vol. 66, pp. 397-413, 1993. 128. R. Pini, M. C. Cavallini, F. Bencini, g. Silverstrini, E. Tonon, W. De Alfieri, N. Marchionni, M. Di Bari, R. B. Devereux, G. Masotti and M. J. Roman, “Cardiovascular remodeling is greater in isolated systolic hypertension than in diastolic hypertension in older adults: the insufficienza cardiaca negli anziani residenti (ICARE) a Dicomano study,” J Am Coll Cardiol, vol. 40, pp. 1283-1289, 2002. 114 References 129. G. L. Pressman, P. M. Newgard, “A transducer for the continuous external measurement of arterial blood pressure” IEEE Trans Biomed Electron, vol. BME-10, pp. 73-81, 1963. 130. M.W. Rajzer, W. Wojciechowska, M. Klocek, I. Palka, M. Brzozowska- Kiszka and K. Kawecka-Jaszcz, “Comparison of aortic pulse wave ve- locity measured by three techniques: Complior, SphygmoCor and Ar- teriograph,” J Hypertens, vol. 26, pp. 2001-2007, 2008. 131. C. H. Reinsch, “Smoothing by spline functions,” Numer Math, vol. 10, 177-183, 1967. 132. G. Rödig, C. Prasser, C. Keyl, A. Liebold and J. Hobbhahn, “Con- tinuous cardiac output measurement: pulse contour analysis vs ther- modilution technique in cardiac surgical patients,” Br J Anaesth, vol. 82, pp. 525-530, 1999. 133. W. J. Rogers, Y.-L. Hu, D. Coast, D. A. Vido, C. M. Kramer, R. E. Pyeritz and N. Reichek, “Age-associated changes in regional aortic pulse wave velocity,” J Am Coll Cardiol, vol. 38, pp. 1123-1129, 2001. 134. M. J. Roman, P. S. Saba, R. Pini, M. Spitzer, T. G. Pickering, S. Rosen, M. H. Alderman and R. B. Devereux, “Parallel cardiac and vascular adaptation in hypertension,” Circulation, vol. 86, pp. 1909- 1918, 1992. 135. M. E. Safar, B. P. Pannier, P. J. Lacolley and B. I Levy, “Cardiac mass and aortic distensibility following calcium blockade in hypertension,” J Cardiovasc Pharmacol, vol. 17(Suppl. 2), pp. S75-S80, 1991. 136. M. E. Safar, J. Blacher, B. Pannier, A. P Guerin, S. J. Marchais, P.-M. Guyonvarch’h, G. M. London, “Central pulse pressure and mortality in end-stage renal disease,” Hypertension, vol. 39, pp. 735-738, 2002. 137. M. E. Safar, B. I. Levy, H. Struijker-Boudier, “Current perspectives on arterial stiffness and pulse pressure in hypertension and cardiovascular diseases,” Circulation, vol. 107, pp. 2864-2869, 2003. 138. M. E. Safar, J. Blacher, A. Protogerou and A. Achimastos, “Arterial stiffness and central hemodynamics in treated hypertensive subjects according to brachial blood pressure classification,” J Hypertens, vol. 26, pp. 130-137, 2008. 115 References 139. K. J. Salaymeh and A. Banerjee, “Evaluation of arterial stiffness in children with William syndrome: does it play a role in evolving hy- pertension?” Am Heart J, vol. 142, pp. 549–555, 2001. 140. G. G. S. Sandor, T. Hishitani, R. E. Petty, M. T. Potts, A. DeSouza, E. DeSouza and J. E. Potts, “A novel Doppler echocardiographic method of measuring the biophysical properties of the aorta in pediatric pa- tients,” J Am Soc Echocardiogr, vol. 16, pp. 745-750, 2003. 141. P. Segers, J. De Backer, D. Devos, S. I. Rabben, T. C. Gillebert, L. M. Van Bortel, J. De Sutter, A. De Paepe and P. R. Verdonck, “Aortic reflection coefficients and their association with global indexes of wave reflection in healthy controls and patients with Marfan’s syndrome,” Am J Physiol Heart Circ Physiol vol. 290, pp. 2385-2392, 2006. 142. P. Segers, E. R. Rietzchel, M. L. De Buyzere, D. De Bacquer, L. M. Van Bortel, G. De Backer, T. C. Gillebert and P. R. Verdonck, “Assessment of pressure wave reflection: getting the timing right,” Physiol Meas vol. 28, pp. 1045-1056, 2007. 143. H. Senzaki, M. Akagi, T. Hishi, A. Ishizawa, M. Yanagisawa, S. Ma- sutani, T. Kobayashi and S. Awa, “Age-associated changes in arterial elastic properties in children,” Eur J Pediatr, vol.161, pp. 547–551, 2002. 144. G. F. Sholler, C. M. Whight, and J. M. Celermajer, “Pulsed Doppler echocardiographic assessment, including use of aortic leaflet separa- tion, of cardiac output in children with structural heart disease,” Am J Cardiol, vol. 57, pp. 1195–1197, 1986. 145. G. F. Sholler, J. M. Celermajer, C. M. Whight, and A. E. Bauman, “Echo Doppler assessment of cardiac output and its relation to growth in normal infants,” Am J Cardiol, vol. 60, pp. 1112–1116, 1987. 146. C. Slagt, J. Beute, M. Hoeksema, I. Malagon, J.-W. R. Mulder and J. A. B. Groeneveld, “Cardiac output derived from arterial pressure waveform analysis without calibration vs. thermodilution in septic shock: evolving accuracy of software versions,” Eur J Anaesthesiol, vol. 27, pp. 550-554, 2010. 147. H. Smulyan, D. S. Siddiqui, R. J. Carlson, G. M. London and M. E. Safar, “Clinical utility of aortic pulses and pressures calculated from 116 References applanated radial-artery pulses,” Hypertension, vol. 42, pp. 150-155, 2003. 148. H. Smulyan, P. R. Sheehe and M. E. Safar ME, “A preliminary evalu- ation of the mean arterial pressure as measured by cuff oscillometry,” Am J Hypertens , vol. 21, pp. 166-171, 2008. 149. S. Söderström, G. Nyberg, M. F. O’Rourke, J. Sellgren and J. Pontén, “Can a clinically useful aortic pressure wave be derived from a radial pressumre wave?” Br J Anaesth, vol. 88, pp. 481-488, 2002. 150. R. W. Stadler, W. C. Karl and R. S. Lees, “New methods for arterial diameter measurement from B-mode images,” Ultrasound Med Biol, vo. 22, pp. 25-34, 1996. 151. R. W. Stadler, J. A. Taylor and R. S. Lees, “Comparison of B-mode, M-mode and echo-tracking methods for measurement of the arterial distension waveform,” Ultrasound Med Biol, vol. 23, pp. 879-887, 1997. 152. C. W. Stetz, R. G. Miller, G. E. Kelly, and T. A. Raffin, “Reliability of the thermodilution method in the determination of cardiac output in clinical practice,” Am Rev Respir Dis, vol. 126, pp. 1001–1004, 1982. 153. P. Studinger, Z. LéNárd, R. Reneman and M. Kollai, “Measurement of aortic arch distension wave with the echo-track technique,” Ultrasound Med Biol, vol. 26, pp. 1285-1291, 2000. 154. K. Sutton-Tyrrell, S. S. Najjar, R. M. Boudreau, L. Venkitachalam, V. Kupelian, E. M. Simonsick, R. Havlik, E. G. Lakatta, H. Spurgeon, S. Kritchevsky, M. Pahor, D. Bauer and A. Newman, “Elevated aortic pulse wave velocity, a marker of arterial stiffness, predicts cardiovas- cuar events in well-functioning older adults,” Circulation, vol. 111, pp. 3384-3390, 2005. 155. G. Swamy and R. Mukkamala, “Estimation of the aortic pressure waveform and beat-to-beat relative cardiac output changes from mul- tiple peripheral artery pressure waveforms,” IEEE Trans Biomed Eng, vol. 55, pp. 1521–1529, 2008. 156. K. Tachibana, H. Imanaka, M. Takeuchi, Y. Takauchi, H. Miyano and M. Nishimura, “Noninvasive cardiac output measurement using par- tial carbon dioxide rebreathing is less accurate at settings of reduced 117 References minute ventilation and when spontaneous breathing is present,” Anes- thesiology, vol. 98, pp.830-837, 2003. 157. H. L. Tan, M. Pinder, R. Parsons, B. Roberts and P. V. van Heerden, “Clinical evaluation of USCOM ultrasonic cardiac output monitor in cardiac surgical patients in intensive care unit,” Br J Anaesth, vol. 94, pp. 287-291, 2005. 158. G. A. Tannenbaum, D. Mathews and C. Weissman, “Pulse contour cardiac output in surgical intensive care unit patients,” J Clin Anesth, vol. 5, pp. 471-478, 1993. 159. B. Trachet, P. Reymond, J. Kips, A. Swillens, M. De Buyzere, B. Suys, N. Stergiopulos and P. Segers, “Numerical validation of a new method to assess aortic pulse wave velocity from a single recording of a brachial artery waveform with an occluding cuff,” Ann Biomed Eng, vol. 38, pp. 876-888, 2010. 160. S. M. Toy , J. Melbin and A. Noordergraaf, “Reduced models of ar- terial systems,” IEEE Trans Biomed Eng, vol. 32, pp. 174-176, 1985. 161. L. M. Van Bortel, E. J. Balkestein, J. J. van der Heijden-Spek, F. H. Vanmolkot, J. A. Staessen, J. A. Kragten, J. W. Vredeveld, M. E. Safar, H. A. Struijker Boudier and A. P. Hoeks, “Non-invasive as- sessment of local arterial pulse pressure: comparison of applanation tonometry and echo-tracking,” J Hypertens, vol. 19, pp. 1037-1044, 2001. 162. J. J. van der Heijden-Spek, J. A. Staessen, R. H. Fagard, A. P. Hoeks, H. A. Struijker Boudier and L. M. Van Bortel, “Effect of age on brachial artery wall properties differs from the aorta and is gender dependent: a population study,” Hypertension, vol. 35, pp. 637-642, 2000. 163. S. Van Huffel and J. Vandewalle, The total least squares problem: com- putation aspects and analysis, Philadelphia: SIAM, 1991. 164. L. E. M. van Lelyveld-Haas, A. R. H. van Zanten, G. F. Borm and D. H. T. Tjan, ”Clinical validation of the non-invasive cardiac output monitor USCOM-IA in critically ill patients,” Eur J Anaesthesiol, vol. 25, pp. 917-924, 2008. 118 References 165. P. Van Overschee and B. De Moor, “N4SID: Subspace algorithms for the identification of combined deterministic-stochastic systems,” Au- tomatica, vol.30, pp. 75–93, 1994. 166. S. J. Vermeersch, E. R. Rietzschel, M. L. De Buyzere, D. De Bacquer, G. De Backer, L. M. Van Bortel, T. C. Gillebert, P. R. Verdonck and P. Segers, “Determining carotid artery pressure from scaled diame- ter waveforms: comparison and validation of calibration techniques in 2026 subjects,” Physiol Meas, vol. 29, pp. 1267-1280, 2008. 167. F. Verbeke, W. Van Biesen, P. Peeters, L. M. Van Bortel and R. C. Vanholder, “Arterial stiffness and wave reflections in renal transplant recipients,” Nephrol Dial Trasplant, vol. 22, pp. 3021-3027, 2007. 168. T. K. Waddell, A. M. Dart, T. L. Medley, J. D. Cameron, B. A. Kingwell, “Carotid pressure is a better predictor of coronary artery disease severity than brachial pressure,” Hypertension, vol. 38, pp. 927-931, 2001. 169. F. J. Walther, B. Siassi, N. A. Ramadan, A. K. Ananda, and P. Y. K. Wu, “Pulsed Doppler determinations of cardiac output in neonates: normal standards for clinical use”, Pediatrics, vol. 76, pp. 829–833, 1985. 170. H. R. Warner, H. J. C. Swan, D. C. Connolly, R. G. Tompkins and E. H. Wood, “Quantitation of beat-to-beat changes in stroke volume from the aortic pulse contour in man,” J Appl Physiol, vol. 5, pp. 495-507, 1953. 171. H. Watanabe, M. Kawai, T. Sibtata, M. Hara, H. Furuhata and S. Mochizuki, “Noninvasive measurement of aortic pressure waveform by ultrasound,” Heart Vessels, vol. 13, pp. 79-86, 1998. 172. K. H. Wesseling, B. de Wit, J. A. P. Weber and N. T. Smith, “A simple device for the continuous measurement of the cardiac output,” Adv Cardiovasc Phys, vol. 5, pp. 16-52, 1983. 173. K. H. Wesseling, J. R. C. Jansen, J. J. Settels, and J. J. Schreuder, “Computation of aortic flow from pressure in humans using a nonlin- ear, three-element model,” J Appl Physiol, vol. 74, pp. 2566–2573, 1993. 119 References 174. W. White, A. Berson, C. Robbins, M. Jamieson, L. Prisant, E. Roc- cella E and S. Sheps, “National standard for measurement of resting and ambulatory blood pressures with automated sphygmomanome- ters,” Hypertension, vol. 21, pp. 504-509, 1993. 175. L. Xie and L. Ljung, “Estimate physical parameters by black-box mod- eling,” in Proc. 21st Chinese Control Conf., 2002, pp. 673-677. 176. K. Yamashita, T. Nishiyama, T. Yokoyama, H. Abe and M. Manabe, “Effects of vasodilation on cardiac output measured by pulseCO,” J Clin Monit Comput, vol. 21, pp. 335-339, 2007. 177. M. Zamir, The Physics of Pulsatile Flow. New York: Springer-Verlag, 2005. 120 Appendix A The Analogy between a Transmission Line and an Elastic Tube with Incompressible Flow Pulsatile flow inside an elastic tube is governed by the Navier-Stokes equa- tions. To derive the most general form of the Navier-stokes equations, two assumptions are made. First, it is assumed that the fluid is Newtonian, which means the velocity gradient inside the fluid is linearly related to the shear stress. It has been a long-standing debate whether blood is a Newto- nian fluid. However, it has been shown that treating blood as a Newtonian fluid is adequate for the study of global dynamics of the circulation and the heart [177]. Second, it is assumed that the flow is laminar, which means the fluid elements only move in the direction of the flow, which is parallel to the axis of the tube. Blood flow is not always laminar inside the cardiovascular circulation. However, at the higher levels of the vascular tree, blood flow is mostly laminar [177]. The derivation of the most general form of the Navier-Stokes equations based on cylindrical coordinates is described in [177]. These equations can be further simplified into the one-dimensional Navier- Stokes equations given below. The derivation is also described in [177]. The equations describing a transmission line are given by: − ∂V (x, t) ∂x = RI(x, t) + L ∂I(x, t) ∂t (A.1) −∂I(x, t) ∂x = C ∂V (x, t) ∂t (A.2) where V is voltage, I is current, R is the line resistance, L is the line induc- tance, C is the line capacitance, and x and t are the spatial and temporal variables respectively. 121 Appendix A. The Analogy between a Transmission Line and an Elastic Tube with Incompressible Flow Equations (A.1) and (A.2) are analogous to the one-dimensional Navier- Stokes equations describing pulsatile incompressible flow in an elastic tube, which are given by [177]: − ∂P (x, t) ∂x = R′Q(x, t) + ρ A ∂Q(x, t) ∂t (A.3) −∂Q(x, t) ∂x = ∂A ∂P ∂P (x, t) ∂t (A.4) where P is pressure, Q is flow, ρ is liquid density, A is the cross-sectional area of the tube, and R′ is a resistive term dependent on the lumen area, and wall and liquid viscosity. By comparing these two sets of equations, it can be seen that voltage corresponds to pressure, current corresponds to flow, line resistance corre- sponds to the resistive term R′, line inductance corresponds to ρ/A, and line capacitance corresponds to ∂A/∂P . In addition, if the transmission line is terminated by a complex load (Figure 2.1), then the terminating condition is given by: ( CL + CL RL RT ) dV (l, t) dt + V (l, t) RP = I + CLRL dI(l, t) dt (A.5) where CL is the load capacitance, RL is the load resistance, RT is the ter- minating resistance and l is the length of the line. 122 Appendix B Implementation of the Transmission Line Model using the Finite Element Method The transmission line equations can be solved numerically using a finite element method [96]. The line is divided into a finite number of elements. Each element has points called interpolation nodes. The voltage and current are evaluated only at a finite number of nodes. Each voltage node is placed between two current nodes and each current node is placed between two voltage nodes. The leap-frog scheme is implemented to evaluate the voltage and current at certain time steps. The voltage is evaluated at half a time step before and after the current. Applying the finite element method to transmission line equations (A.1) and (A.2) gives the following equations: I n+ 1 2 j = K1I n− 1 2 j −K2(V ni+1 − V ni ) (B.1) V n+1j = K3V n i −K4(I n+ 1 2 j+1 − I n+ 1 2 j ) (B.2) where i and j indicate the nodes and n indicates the time step. K1 to K4 are given by: K1 = 1 (B.3) K2 = ∆t L∆x (B.4) K3 = 1 (B.5) K4 = ∆t C∆x (B.6) For numerical stability, we require u ≤ ∆x ∆t (B.7) 123 Appendix B. Implementation of the Transmission Line Model using the Finite Element Method where u is the wave propagation speed. The terminating condition (A.5) is approximated as follows: A5V n+1 m = A3I n+ 1 2 p +A4 I n+ 1 2 p − In− 1 2 p ∆t +A6V n m (B.8) where m indicates the last voltage node and p indicates the last current node. The last current node is located before the last voltage node. A1 to A6 are given as follows: A1 = CL +CL RL RT (B.9) A2 = 1 RT (B.10) A3 = 1 (B.11) A4 = CLRL (B.12) A5 = A1 ∆t + A2 2 (B.13) A6 = A1 ∆t − A2 2 (B.14) 124 Appendix C The Derivation of the Transfer Function from the Ascending Aortic Flow to the Radial Artery Pressure For a uniform lossless transmission line, Po = Pi + Pr (C.1) Qo = 1 Zch (Pi − Pr) (C.2) Pt = Piz −d + Prz d (C.3) where Po is the pressure at the beginning of the line, Pi is the incident wave at the beginning of the line, Pr is the reflected wave at the beginning of the line, Qo is the flow at the beginning of the line, Zch is the characteristic impedance of the line, Pt is the pressure at the end of the line, n is the discrete time delay for the wave to travel from the beginning to the end of the line, z is the z-variable. Pi and Pr are related by the wave reflection coefficient at the tube-load interface Γ as: Γ = Prz d Piz−d (C.4) Using equations (C.2)-(C.4), the transfer function (M(z)) from Qo to Pt is then given by: M(z) = Pt Qo = Zch(1 + Γ)z n z2n − Γ (C.5) Γ is given by: Γ = ZL − Zch ZL + Zch (C.6) 125 Appendix C. The Derivation of the Transfer Function from the Ascending Aortic Flow to the Radial Artery Pressure where ZL is the impedance of the terminating load. In terms on the s-variable (continuous-time), ZL is given by: ZL = 1 + CLRLs 1 RT + (CL + CLRL RT )s (C.7) The s-variable is replaced by the z-variable (discrete-time) using the follow- ing approximation: s = 2Fs(z − 1) (z + 1) (C.8) where Fs is the sampling frequency. Then ZL is given by: ZL = AZz +BZ CZz +DZ (C.9) where AZ = 1 + 2FsCLRL (C.10) BZ = 1− 2FsCLRL (C.11) CZ = 1 RT + 2Fs(CL + CL RL RT ) (C.12) DZ = 1 RT − 2Fs(CL + CLRL RT ) (C.13) Using equations (C.6), (C.9)-(C.13), Γ can be written as: Γ = Arz +Br Crz +Dr (C.14) where Ar = AZ − ZchCz (C.15) Br = BZ − ZchDz (C.16) Cr = AZ + ZchCz (C.17) Dr = BZ + ZchDz (C.18) Using equations (C.5), (C.14)-(C.18), M(z) can be written as: M(z) = AMz d+1 +BMz d CMz2d+1 +DMz2d + EMz + FM (C.19) 126 Appendix C. The Derivation of the Transfer Function from the Ascending Aortic Flow to the Radial Artery Pressure where AM = Zch(Ar + Cr) (C.20) BM = Zch(Br +Dr) (C.21) CM = Cr (C.22) DM = Dr (C.23) EM = −Ar (C.24) FM = −Br (C.25) In terms of CL, RL and RT , equations (C.20)-(C.25) are given by: AM = 2Zch(1 + 2FsCLRL) (C.26) BM = 2Zch(1− 2FsCLRL) (C.27) CM = Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.28) DM = Zch RT + 1− 2FsCL(RL + Zch + RLZch RT ) (C.29) EM = Zch RT − 1− 2FsCL(RL − Zch − RLZch RT ) (C.30) FM = Zch RT − 1 + 2FsCL(RL − Zch − RLZch RT ) (C.31) M(z) can also be written in the more common form as M(z) = b1z −d + b2z −d−1 1 + a1z−1 + a −2d 2 + a −2d−1 3 (C.32) with b1 = 2Zch(1 + 2FsCLRL) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.33) b2 = 2Zch(1− 2FsCLRL) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.34) a1 = Zch RT + 1− 2FsCL(RL + Zch + RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.35) a2 = Zch RT − 1− 2FsCL(RL − Zch − RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.36) a3 = Zch RT − 1 + 2FsCL(RL − Zch − RLZchRT ) Zch RT + 1 + 2FsCL(RL + Zch + RLZch RT ) (C.37) 127"@en ; edm:hasType "Thesis/Dissertation"@en ; vivo:dateIssued "2011-11"@en ; edm:isShownAt "10.14288/1.0103216"@en ; dcterms:language "eng"@en ; ns0:degreeDiscipline "Electrical and Computer Engineering"@en ; edm:provider "Vancouver : University of British Columbia Library"@en ; dcterms:publisher "University of British Columbia"@en ; dcterms:rights "Attribution-NonCommercial-NoDerivatives 4.0 International"@en ; ns0:rightsURI "http://creativecommons.org/licenses/by-nc-nd/4.0/"@en ; ns0:scholarLevel "Graduate"@en ; dcterms:title "The measurement and the monitoring of the central hemodynamics using non-invasive and minimally-invasive data"@en ; dcterms:type "Text"@en ; ns0:identifierURI "http://hdl.handle.net/2429/35847"@en .