## 1 Introduction  
Numerical models have long been employed to better understand the complex physics involved in human phonation. For example, reduced-order and finite element numerical models of the vocal folds (VFs) can self-oscillate in a manner representative of actual VF kinematics during sustained vowels [1], pitch glides [2], and, in a few cases, running speech [3]. Such models have explored a wide range of phenomena relevant to normal and pathological phonation, including the impact of a posterior glottal gap [4], the ventricular folds [5], phonation onset pressure [6], and the efficacy of various compensation mechanisms for vocal hyperfunction [7; 8]. Whereas the majority of modeling efforts in speech, to date, have employed models with general population-based parameters to uncover the universal physical underpinnings of human phonation, a few research teams have begun to explore the development of subject-specific numerical models for phonation [9; 10; 11; 12; 13]. The dynamics of the VFs are sensitive to a variety of factors, including subglottalpressure [14], laryngeal muscle activation [15], and a posterior glottal gap [4], to name a few. These factors, which can play a role in vocal hyperfunction [16] and other pathologies [17; 18; 19; 20], can be a challenge to observe clinically. Subject-specific models, on the other hand, are constructed based on measurements of the subject through less challenging media, such as high speed videoundoscopy (HSV) [9; 11; 21] and offer the potential to elucidate clinically opaque features and parameters, such as VF contact pressures. In general, development of a subject-specific model entails estimation of numerical parameters such that the model behavior \"matches\" some observed data from the individual of interest, such as VF kinematics from HSV [9; 11; 21]. Specifically, numerical model parameters are sought such that some model output(s) (e.g., the VF kinematics) best match the equivalent observation data from the subject. This problem is ill-posed, thus necessitating inverse analysis techniques to determine the model parameters [22]. To date, optimization [9; 22], machine learning [23], and Bayesian [10; 11; 21] frameworks have been developed to estimate VF model parameters. Optimization-based approaches define a cost functional and employ optimization techniques to determine the parameter values that minimize the functional. [PERSON] et al. [22] was the first to successfully use this approach, employing the Nelder-Mead algorithm to determine the vibrating masses, spring stiffnesses, and subglottal pressure of a two mass VF model by minimizing the least-squares error between specific Fourier coefficients of the measured VF trajectories and the simulated trajectories. Genetic algorithms were subsequently employed to determine parameters of a two mass model that best reproduced the trajectories from patients suffering from unilateral VF paralysis [16]. Inverse procedures have also been used to classify disordered versus healthy VF oscillations. [PERSON] et al. [24] employed simulated annealing to minimize the Euclidean distance between the model and experimental vocal fold trajectories to distinguish between normal and dysphonic subjects during sustained vowel phonation and a pitch glide. Very recently, deep learning tools have been employed for subject-specific modeling [23]. A long short-term memory network, trained on simulated data, was used to estimate the subglottal pressure from ex vivo recordings of porcine vocal folds. They were able to produce accurate estimates of the subglottal pressure with very low online computational costs. The vast majority of these subject-specific VF modeling efforts, and all of those mentioned above, have used lumped-element representations of the VFs and employed inverse methods to determine the reduced-order parameters. Of noted exception are the efforts by [PERSON] et al. [12] and [PERSON] et al. [13] that employ finite element model representations of the VFs. [PERSON] et al.