## 2 Materials and Methods
### Development of the Musculoskeletal Model of Swallowing
#### 2.2.5 Muscle Force Estimation Using an Optimization Program  
(1) Equation of Motion  
Using Lagrangian mechanics and setting the generalized coordinates as \(\mathbf{q}\), the dynamic equations (such as Equation (8)) of the musculoskeletal model in the following form were obtained:  
\[\mathbf{M\ddot{q}}+\mathbf{C\dot{q}}+\mathbf{Kq}\ +\ \mathbf{F}_{pu}^{m}+\mathbf{g}(\mathbf{q})\ =\ \mathbf{Q} \tag{8}\]  
where \(\mathbf{M}\) is the inertia matrix with a mass \(M\) and a moment of inertia \(J\) and Cand Kare damping and stiffness matrices, respectively. The \(\mathbf{F}_{pu}^{m}\) is the matrix representing the passive force of the muscles. The matrix \(\mathbf{g}\) represents the effect of gravity. The \(\mathbf{q}=\left[x^{h},\ y^{h},\ z^{h},\ \theta^{h},\ x^{th},\ y^{h},\ z^{th},\ \theta^{th}\right]^{T}\) stands for positions and postures of the hyoid bone and thyroid cartilage. The \(\mathbf{Q}=\left[F_{x}^{h},\ F_{y}^{h},\ F_{z}^{h},\ T^{h},\ F_{x}^{th},\ F_{y}^{ th},\ F_{z}^{th},\ T^{th}\right]^{T}\) is the generalized force, which includes forces and torques acting on the hyoid bone and thyroid cartilage. The \(h\) and \(th\) indices stand for parameters related to the hyoid bone and thyroid cartilage, respectively. The mass of the hyoid bone and thyroid cartilage, \(M^{h}\) and \(M^{th}\), were estimated to be \(5.0\times 10^{-3}\) kg and \(1.5\times 10^{-2}\) kg, respectively, by calculating the volumes from the skeleton data. Further, the moment of inertia about the center of gravity of the hyoid bone and thyroid cartilage, \(J^{th}\) and \(J^{th}\), were estimated using the function of 3D CAD software, SoldWorks(r), to be \(4.12\times 10^{-7}\) kg m\({}^{2}\) and \(2.92\times 10^{-6}\) kg m\({}^{2}\), respectively. Using the dynamic Equation (8), inverse dynamics analysis was performed to compute the generalized force \(\mathbf{Q}\) incorporating the forces and torque that realize the motion of the hyoid bone and thyroid cartilage, as shown in Figure 2. Here, the differentiation of \(\mathbf{q}\) with time was computed by a finite-difference method. In addition, the effect of gravity was eliminated (\(\mathbf{g}\ =\ \mathbf{0}\)) for simplicity. (2) Muscle Force Estimation  
The developed musculoskeletal model has a redundant number of muscles compared with the degrees of freedom of the hyoid bone and thyroid cartilage. Due to this redundancy, estimating the muscle forces from the obtained generalized forces was an ill-posed problem, and to solve it, Crowninshield's cost function \(O(f)\)[16] was used. The function was formed as a quadratic program as below:  
\[O(f)=\sum_{i=1}^{N}\left(\frac{f_{i}}{PCSA_{i}}\right)^{2}\to min \tag{9}\]  
\[\text{\emph{subject to}}:\ F_{x}=\sum_{i=1}^{N}f_{ix},\ \ F_{y}=\sum_{i=1}^{N}f_{iy},\ \ T_{\theta}=\sum_{i=1}^{N}f_{i}r_{i} \tag{10}\]  
\[f_{i}^{min}\leq\ f_{i}\ \leq\ f_{i}^{max}\ \ \ \ \ (i\ =\ 1,\cdots,N) \tag{11}\] \begin{table}
\begin{tabular}{c c c} \hline
**Parameter** & **Value** & **Unit** \\ \hline \(\sigma_{max}\) & \(0.7\times 10^{6}\) & [N/m\({}^{2}\)] \\ \(l_{m\text{s}0}\) & Derived from the model & [m] \\ \(l_{c0}\) & \(l_{mus0}\) & [m] \\ \(l_{csh}\) & \(0.5l_{c0}\) & [m] \\ \(V_{vm}\) & \(6l_{c0}\) & [m/s] \\ \(V_{sh}\) & \(0.3\) & [-] \\ \(V_{shl}\) & \(0.23\) & [-] \\ \(V_{ml}\) & \(1.3\) & [-] \\ \(V_{cr}\) & \(0.2\) & [-] \\ \(P_{E\text{\emph{s}}l}\) & \(4.0\) & [-] \\ \(P_{E\text{\emph{s}}m}\) & \(0.4\) & [m] \\ \hline \end{tabular}
\end{table} 
Table 2: Parameters of the muscle model. where \(PCSA_{i}\) is the physiological cross-sectional area (PCSA) of the muscle \(i,f_{i}\) is the muscle force, \(r_{i}\) is the moment arm of muscle \(i\) at the hyoid bone or thyroid cartilage, \(f_{i}^{max}\) is the maximum force of the muscle \(i\), and \(f_{i}^{min}\) is the minimum. A combination of possible muscle forces that minimize the total muscle force was obtained by optimizing the cost function (9). Function (10) is the equality constraint and possible muscle forces were required to satisfy the forces and the torque computed by inverse dynamics.