## Appendix A  
In this appendix, the sensitivity analysis of the state variables of model (1) with respect to the model parameters will be carried out. Consider the general nonlinear \(n\)-th order model defined by  
\[\frac{\mathrm{d}x(t,\mathbf{\theta})}{\mathrm{d}t}=f(x,\mathbf{\theta})+Bu(t), \tag{11}\]  
where \(x(t,\mathbf{\theta})\in\mathbb{R}^{n}\), \(f(x,\mathbf{\theta}):\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}\), \(B\in\mathbb{R}^{n}\), \(u(t)\in\mathbb{R}\) and \(\mathbf{\theta}\in\mathbb{R}^{p}\) is the parameter vector. Now consider function \(S_{ij}(t,\mathbf{\theta}):\mathbb{R}\rightarrow\mathbb{R}\)  
\[S_{ij}(t,\mathbf{\theta})=\frac{\mathrm{d}x_{i}(t,\mathbf{\theta})}{\mathrm{d}\theta_{ j}}, \tag{12}\]  
which defines the sensitivity of the \(i\)-th state variable \(x_{i}(t,\mathbf{\theta})\) to the variations in the \(j\)-th parameter \(\theta_{j}\). The time derivative \(\frac{\mathrm{d}S_{ij}(t,\mathbf{\theta})}{\mathrm{d}t}\) of the sensitivity function \(S_{ij}(t,\mathbf{\theta})\) can be derived as  
\[\frac{\mathrm{d}S_{ij}(t,\mathbf{\theta})}{\mathrm{d}t}=\frac{\mathrm{d}f_{i}(x,\bm {\theta})}{\mathrm{d}\theta_{j}}+\sum_{k=1}^{n}A_{ik}(x,\mathbf{\theta})S_{kj}(t, \mathbf{\theta}), \tag{11}\]  
where \(A_{ik}(x,\mathbf{\theta})\) are the entries of the Jacobian matrix (37). By numerically solving differential equations \(\frac{\mathrm{d}S_{ij}(t,\mathbf{\theta})}{\mathrm{d}t}\) for \(i=1,2,\ldots,n\), and \(j=1,2,\ldots,p\) assuming \(S_{ij}(0,\mathbf{\theta})=\mathbf{0}\) simultaneously with the differential equations of model (11), the sensitivity functions \(S_{ij}(t,\mathbf{\theta})\) can be found. Considering \(f(x,\mathbf{\theta})\) defined by (1) and \(\mathbf{\theta}\) defined by (3), we can find the nonzero partial derivatives \(\frac{\mathrm{d}f_{i}(x,\mathbf{\theta})}{\mathrm{d}\theta_{i}}\) as  
\[ \begin{array}{c}\frac{\mathrm{d}f_{1}(x,\mathbf{\theta})}{\mathrm{d}\theta_{1}} =N(t)(1-bN(t)),\frac{\mathrm{d}f_{1}(x,\mathbf{\theta})}{\mathrm{d}\theta_{1}}=- aN^{2}(t),\,\frac{\mathrm{d}f_{3}(x,\mathbf{\theta})}{\mathrm{d}\theta_{3}}=T(t)(1- dT(t)),\,\frac{\mathrm{d}f_{3}(x,\mathbf{\theta})}{\mathrm{d}\theta_{4}}=-cT^{2}(t),\\ \frac{\mathrm{d}f_{2}(x,\mathbf{\theta})}{\mathrm{d}\theta_{5}}=T(t)N(t),\,\frac{ \mathrm{d}f_{2}(x,\mathbf{\theta})}{\mathrm{d}\theta_{6}}=-L(t),\,\frac{\mathrm{d }f_{1}(x,\mathbf{\theta})}{\mathrm{d}\theta_{7}}=-T(t)N(t),\,\frac{\mathrm{d}f_{3 }(x,\mathbf{\theta})}{\mathrm{d}\theta_{8}}=-T(t)N(t),\\ \frac{\mathrm{d}f_{2}(x,\mathbf{\theta})}{\mathrm{d}\theta_{9}}=-L(t)T(t),\frac{ \mathrm{d}f_{3}(x,\mathbf{\theta})}{\mathrm{d}\theta_{10}}=-u(t),\,\frac{\mathrm{d }f_{1}(x,\mathbf{\theta})}{\mathrm{d}\theta_{12}}=-N(t)u(t),\\ \frac{\mathrm{d}f_{2}(x,\mathbf{\theta})}{\mathrm{d}\theta_{13}}=-L(t)u(t),\,\frac{ \mathrm{d}f_{3}(x,\mathbf{\theta})}{\mathrm{d}\theta_{14}}=-T(t)u(t).\end{array}  \tag{12}\]  
Now, considering the initial conditions given by (47), the state vector of the nominal model was simulated and \(S_{ij}(t,\mathbf{\theta})\) were determined by solving (11). The evolution of individual sensitivities \(S_{ij}(t,\mathbf{\theta})\) scaled with respect to the given nominal parameter value \(\overline{\theta}_{j}\) to obtain the sensitivity \(S_{ij}(t,\mathbf{\theta})\overline{\theta}_{j}\) to relative variation in a parameter \(\frac{\theta_{j}}{\overline{\theta}_{j}}\) is visualized in Figures 11-13. By analyzing the sensitivities \(S_{ij}(t,\mathbf{\theta})\) visualized in Figures 11-13, one can make the following conclusions concerning the identifiability of the individual parameters. Sensitivities of all three state variables \(N(t)\), \(L(t)\), \(T(t)\) to the variations in the parameters \(\theta_{4}\) (denoted \(d\)), \(\theta_{8}\) (denoted \(\alpha_{2}\)), \(\theta_{9}\) (denoted \(\beta_{1}\)) are consistently very low compared to typical magnitudes of the states.