Respiratory Sinus Arrhythmia (RSA) is a change in the heart rate that corresponds to the frequency of respiration, but its causative mechanisms in humans remain only partially determined. In this thesis, models of the human cardiovascular system have been developed to give physiologically reasonable explanations of RSA. In a normally intact system, the cardiovascular system has respiration as input and heart rate as the output, which can be used to evaluate the cardiovascular system dynamics. It is proposed that respiration oscillations play a major role in the generation of RSA in healthy humans through the action of medullary respiratory neurons, which have been shown to both control respiration itself and via a coupling to autonomic centers to modulate heart rate by varying parasympathetic and sympathetic inputs to the sinoatrial node. The Baroreflex does not play an important role in modifying the RSA response in normal conditions but does maintain arterial blood pressure within a fairly narrow range through feedback control. Based on the neural network autoregressive with exogenous input (ARX) model, sensitivity analysis has been employed to demonstrate that instantaneous lung volume (ILV) has a much higher percentage of contribution to heart rate variability (HRV) than systolic blood pressure (SBP) (typically 11.6%, 4.6%), which agrees well with our cardiovascular system model. System identification of RSA in random, regular and spontaneous breathing patterns performed well when using an ARX model of low order, typically [4, 2, x]. In order to get a tradeoff between the accuracy of the model and an excessive number of parameters associated with that model, both the false nearest neighbors (FNN) algorithm and visual inspection of the loss function were used for model order determination. Good model qualities were proved by model validations and agreed well with the results of previous research. Studying transfer characteristics from fluctuations of ILV to HRV, the impulse response was obtained from the ARX model and it was then decomposed algebraically into two combinations of exponential decays, i.e., the fast and slow response components. The low order ARX model was tested on data acquired from a range of healthy subjects and patients. Analysis of the data indicated that the fast response component of the impulse response corresponds to the high frequency power (HF) (0.3-0.5 Hz) of HRV fluctuations, and the slow response component corresponds to the low frequency power (LF) (0.08-0.15 Hz) of HRV fluctuation. The fast component reflects the change of HRV due to parasympathetic input and the slow component is attributed to both sympathetic and parasympathetic input. Analysis of the model parameters, the amplitude and time constant of impulse response components, showed these can identify differences in the system response due to the shifts in autonomic balance produced by changing posture and in patients with large sympathetic inputs. Compared with HF:LF power ratio, our slow components of impulse response are more informative in quantifying the balance between the sympathetic and parasympathetic response due to respiration. To our knowledge, this is first attempt to use an ARX model to quantify RSA for patients in uncontrolled situations where the spontaneous breathing pattern can take on a variety of forms. Our study shows that the results of this system identification to patients are reliable quantitative indexes of RSA after the data has been examined to ensure that cardiorespiratory interaction is linear and stationary.