Biomedical and biophysical limits to mathematical modeling of pulmonary system mechanics: a scoping review on aerosol and drug delivery.

Undoubtedly, the construction of the biomechanical geometry systems with the help of computer tomography (CT) and magnetic resonance imaging (MRI) has made a significant advancement in studying in vitro numerical models as accurately as possible. However, some simplifying assumptions in the computational studies of the respiratory system have caused errors and deviations from the in vivo actual state. The most important of these hypotheses is how to generate volume from the point cloud exported from CT or MRI images, not paying attention to the wall thickness and its effect in computational fluid dynamic method, statistical logic of aerosol trap in software; and most importantly, the viscoelastic effect of respiratory tract wall in living tissue pointed in the fluid-structure interaction method. So that applying the viscoelastic dynamic mesh effect in the form of the moving deforming mesh can be very effective in achieving more appropriate response quality. Also, changing the volume fraction of the pulmonary extracellular matrix constituents leads to changes in elastic modulus (storage modulus) and the viscous modulus (loss modulus) of lung tissue. Therefore, in the biomedical computational methods where the model wall is considered flexible, the viscoelastic properties of the texture must be considered correctly.

Introduction

So far, immense investigations have been done on the respiratory system function and different perspectives of the human breathing procedure. Many researchers investigated the simulation of the airflow and deposition of the particle via the diffusion procedure in a 3D human airway model (Zhang and Lessmann 1996; Cebral and Summers 2004). Also, they presented that both the outlet and inlet boundary conditions of the flow have considerable results in the flow movement. Furthermore, they investigated the aerodynamics of the flow in the human tracheal and bronchi using CFD. Also, raised shear stress and reduced pressure in the stenosis zone were demonstrated. Moreover, the accurate choice of the boundary condition of the outflow in an actual sneeze can considerably affect the stenosis zone conditions as it also analyzes and confirmed in the previous studies (Mortazavy Beni et al. 2019a, b, 2021; Mortazavi et al. 2020). Images produced from CT with formatting DICOM (Digital Imaging and Communications in Medicine) were imported to Mimics software using thresholding to reconstruct the human respiratory point cloud. A respiratory geometry extracted from the point cloud can be seen in Fig. 1 for a healthy 30-year-old male subject. Most prior studies on the respiratory system have been done applying the CFD method, but indeed, the respiratory mechanism should be investigated as the interaction of the air and the tract walls, which is modeled by the FSI technique.

Fig. 1

A respiratory geometry with the numerical grid generation for a healthy 30-year-old male subject. a A respiratory geometry. b Computational grid generation

Part 1: fine solid particles deposition

A wide range of mathematical models analyzed the airflow and respiratory aerosol deposition to the extrathoracic airways (Xiao, et al. 2021; Xiao, et al. 2021; Saha, et al. 2019; Islam et al. 2021a, b). The availability of the high-resolution CT-scan data for the upper airway model helps the researcher to perform the mathematical analysis of the particle deposition at the extrathoracic airways. On the contrary, a limited number of published studies are available for airflow and particle transport to the lower airways. The lack of high-resolution CT images and high computational cost restricted the researchers from analyzing the lower airways' aerosol deposition. Despite the high computational cost, several studies analyzed the monodisperse (Islam et al. 2017a), polydisperse (Islam et al. 2018), and SARS CoV-2 (Islam et al. 2021c) aerosol to the lower airways of a large-scale model. The numerical study (Islam et al. 2017a) employed a highly complex 17-generations model and analyzed the monodisperse particle transport to lower airways. The computational study reports a higher flow (%) distribution to the right lung than the left lung. However, in reality, the inhaled particles are polydisperse, and the polydisperse aerosol study (Xiao, et al. 2021) reports a highly complex deposition pattern for five different lobes of the whole lung model. A precise knowledge of the fine aerosol transport to the lower airways is important as they can transport to the deeper bronchioles due to their tiny size. Recently, a large-scale 17-generation model performed the ultrafine aerosol transport and deposition to the lower airways (Islam et al. 2017b). The numerical study critically analyzed the flow behavior and calculated a comprehensive generation to generation pressure variation, which is important for the health risk assessment of respiratory patients. The computation study reports lobe-specific deposition patterns for different physical activity conditions. A numerical study investigated the lung burden and ultrafine aerosol deposition to the airways (Salma et al. 2015). A wide range of aerosol diameter (6–1000 nm) is tested, and the numerical study reports higher deposition concentration in tracheobronchial airways than the extrathoracic airways. The study also reports that the deposition concentration at the alveolar region increased with physical activity. A numerical study performed the combined nanosize, and microparticle deposition is to 16-generation airways. The single bifurcation airway model consists of a CT-based realistic (first 4-generation) and non-realistic (5–16 generations) model. The computational study predicts 0.2–0.5% nanocarrier deposition for the single bifurcation model. Breathing has been regarded as the most significant process in the body (Sohrabi et al. 2017). Table 1 summarizes some previous studies with their method and study limitation.

Table 1

Fine particles deposition in the respiratory system

Method	Citation	Strong points of the result	Study limitation
CFD	Zhang et al. (2019)	Aerodynamic characteristics and deposition of dust mite allergens in the nasal cavity were analyzed	Nasal cavity with rigid wall, steady state airflow
CFD	Yan et al. (2019)	The morphological variation of pharynx significantly affects the particle deposition features	Pharynx with rigid wall
CFD	Phoung et al. (2020)	Extrapolating the test results is complicated by anatomical and physiological differences between animals and humans	Larynx with rigid wall, steady airflow
CFD	Djupesland (2013)	Estimate the regional deposition in upper airway and concern to the inlet air and fine solid particles	Upper airway with rigid wall
CFD-PIV	Xu et al. (2020)	Flow field characteristics in the trachea region were measured	Trachea with rigid wall, steady state airflow
CFD	Phuong et al. (2018)	The deposition fractions of the monkey's numerical airway model agreed well with in-vitro and human model data when equivalent Stokes numbers were compared, based on the minimum cross-sectional area as representative of length scale	Upper airway with rigid wall
SLA	Valtonen et al. (2020)	The results in vivo were higher than the results in vitro in maxillary sinus volumes with a ratio of 1.05 ± 0.01 (mean ± SD) and in the nasal cavities with a ratio of 1.20 ± 0.1 (mean ± SD)	Nasal cavity with rigid wall
CFD- SLA	Zhan et al. (2019)	In the middle and upper nasal tract, vortex line separation occurs and there is significant passage effect	Upper respiratory tract with rigid wall
SLA	Kelly et al. (2004a,b)	Information on the deposition efficiency of aerosol particles in the nasal airways is used for optimizing the delivery of therapeutic aerosols into the nose for risk assessment of toxic airborne pollutants inhaled through the nose into the respiratory system	Nasal airways with rigid wall
CFD	Regard Rahimi-Gorji et al. (2016)	Enhancing inhalation flow rate and particle size will largely increase the inertial force and consequently, more particle deposition is evident suggesting that inertial impaction is the dominant deposition mechanism in tracheobronchial airways	Tracheobronchial airways with rigid wall, steady state airflow
CFD- SLA	Collier et al. (2018)	A turbulent laryngeal jet flow was observed and affected remarkably the velocity profiles in the trachea	Laryngeal airway with rigid wall, steady state airflow
CFD	Lieber and Zhao (1998)	The results suggest that under the conditions studied a quasisteady flow assumption for oscillatory flow is valid for only about 50% of the oscillatory period, or it is limited to represent the oscillatory flow only in the vicinity of peak inspiration and peak expiration	Airway tract with rigid wall
CFD	Zhang et al. (2002)	The deposition rates of the steady state analysis are much lower in comparison to their unsteady counterparts	Upper airway with rigid wall
CFD	Bahmanzadeh et al. (2016)	For breathing under a rest condition with a frequency of 0.25  Hz, the quasi-steady airflow assumption in the nasal cavity was found to be reasonable when the instantaneous Strouhal number was smaller than 0.2	Nasal cavity with rigid wall
CFD	Cui et al. (2020)	The properties of airflow structures are highly impacted by the respiration pattern	Upper airway with rigid wall
CFD	Gu et al. (2019)	The total deposition of micro particles ranging from 1 to 20 μm under unsteady inhalation was almost the same as that at steady state when the volume of inhaled airflow was equivalent	Nasal cavity with rigid wall
CFD	Kiasadegh et al. (2020)	Assuming the steady flow could be accurate in predicting the deposition of fibrous particles. Nonetheless, this is not valid for regional particle deposition	Upper airway with rigid wall

Part 2: fine liquid droplets deposition

Liquid droplets transport to the airways attracts the attention of researchers in recent days. It is evident that the SARS CoV-2 virus transmits as droplets during sneezing and coughing (Meselson 2020). A numerical study analyzed the cough-generated SARS CoV-2 virus-laden droplets transport to the airways (April Si et al. 2021). The study used 0.1–4 µm diameter droplets, and the single path model predicts the high exhalation rate for 2-µm diameter droplets. A numerical study investigated the cough-generated droplets' transport and deposition to the airways (Madas, et al. 2020). The computational study predicts significantly higher deposition concentration (seven times) at the upper airways than the acinar region. A recent study performed the SARS CoV-2 aerosol transport to the age-specific airway model (Islam, et al. 2021d). The numerical study calculated the deposition efficiency for 50-, 60- and 70-year-old airways model. The study reports a higher deposition concentration at the right lung for smaller SARS CoV-2 aerosol than the left lung. However, an opposite scenario is observed for larger diameter SARS CoV-2 aerosol. Also, adequate air conditioning in the nasal airways is mandatory for respiration and gas exchange in the lower respiratory tract. A high increase of humidity and temperature at the end of inspiration, in relation to the environmental conditions, was found in the anterior nasal segment (Drettner et al. 1977; Keck et al. 2000). In this regard, research indicates that homeopathic remedies (a) contain measurable source and silica nanoparticles heterogeneously dispersed in colloidal solution; (b) act by modulating the biological function of the allostatic stress response network; (c) evoke biphasic actions on living systems via organism-dependent adaptive and endogenously amplified effects; (d) improve systemic resilience (Mygind and Vesterhauge 1978; Bell and Koithan 2012; Rissler et al. 2012; Wichers et al. 2006). Nasal drug delivery that exploits the olfactory and trigeminal neuronal pathways to deliver drugs to the brain is being widely explored by pharmaceutical companies for the delivery of challenging drugs. Low-molecular-weight and lipophilic drugs are effectively absorbed by the intranasal route for efficacious brain targeting; however, high-molecular-weight and hydrophilic drugs present challenges in intranasal delivery (Bahadur and Pathak 2012; Xi et al. 2015). For drug delivery to the maxillary sinus, particles in the 0.67- to 0.99-μm range had improved efficiency of deposition in the maxillary sinus compared with larger particles after maxillary antrostomy. Larger particles appeared to deposit directly in the nasal vault, while smaller particles were more likely to reach the maxillary sinus (Hilton et al. 2008; Möller et al. 2014). Hence, according to part 1 and part 2 in this section, it can be inferred that the method features can have an impact on aerosol deposition and transport in different regions of the respiratory system. Therefore, in this study, we tried to improve the methods and materials for future perspectives, which is the main motivation of this study.

Improvement of methods and materials

Biomechanical geometry construction with the help of CT and MRI has made a significant advancement in studying in vitro numerical models. However, some simplifying assumptions in the computational studies of the respiratory system have caused errors and deviations from the in vivo actual state. For example, in CFD modeling, the respiratory tract wall is considered rigid, or in FSI modeling, the wall is often considered hyperelastic; although, in reality, the tract wall has a viscoelastic behavior. So that applying the viscoelastic dynamic mesh effect in the form of the MDM can be very effective in achieving a more appropriate response quality. Therefore, considering a viscoelastic airway wall during transient airflow responding to the pressure pulse from the lung has a special attitude. The elastic behaviors of the ECM were simulated using a spring network model and various geometric arrangements (Kowe et al. 1986; Kimmel and Budiansky 1990). However, it is understood that a hexagonal pattern is most appropriate to represent alveolar morphology (Kononov et al. 2001). Therefore, a hexagonal spring arrangement with each spring representing an alveolar wall was utilized. Also, parallel springs provide stability to the network and represent the role of PGs. A spring was placed in parallel with a spring–dashpot to represent the combined effect of PGs, elastin, and collagen (Fig. 2a). Each alveolar wall spring showed the combined mechanical property of collagen, elastin, and the surrounding matrix inside the lung parenchymal cell. Elastin, collagen, and PGs were assumed to behave as a linear spring, and the surrounding matrix was showed a Maxwell element to provided viscoelastic behaviors. According to Fig. 2a, the standard linear solid (SLS) model represents the ECM unit. Also, alveolar wall elements with a viscoelastic model representing elastin, , collagen, , and the surrounding matrix, and . Moreover, PGs elements mimic the role of a linear spring model for the strength of the lung parenchymal cell then . Also, Fig. 2b shown a viscoelastic network representing the generalized viscoelastic model for the normal lung ECM strip ().

Fig. 2

The ECM viscoelastic behavior; a the SLS model for the ECM unit; b the generalized viscoelastic model for the normal lung ECM strip

Viscoelastic properties were included in this spring network employing a spring–dashpot element (Cavalcante et al. 2005; Sopakayang and Vita 2011). In the current study, within a viscous matrix using collagen, elastin, and PGs; hence, a linear spring–dashpot system in parallel with alveolar walls to represent a surrounding matrix. Elastin fibers represent linear elastic property, and the stress within elastin fibers can be determined as ; where is the elastin Young’s modulus and is the axial strain in the model (Sugihara et al. 1971). The employment of collagen in the simulation was modeled using a Weibull cumulative density function (Fung 1967). Also, collagen type I was considered and found the most plentiful and the major tensile type of collagen in ECM (Asgari et al. 2017). Furthermore, the linear approximation model simulated the stress within collagen fibers can be written by the nonlinear equation as ; where is the collagen Young’s modulus inside the collagen fibers (Sugihara et al. 1971). Stabilizing elements exhibiting the role of PGs inside the alveolar walls were contained using a linear spring element with stress determined as ; where is the PGs Young’s modulus (Yuan et al. 2000; An et al. 2007). The stress within the surrounding matrix was represented using a spring–dashpot (Maxwell) model:

Here, and represent the surrounding matrix Young’s modulus and the damper coefficient of the Maxwell model component (Sugihara et al. 1972), respectively. Also, and are the strains within the spring and dashpot, respectively. The total stress across these SLS models (Fig. 2a) was determined according to the rule of a mixture. Therefore, the total stress for the alveolar element can be written as:

And finally, by Eq. 1–3, the constitutive relation for the three elements model shown in Fig. 2a is derived as follows:

The ratio has the units of time and can be remarked as the characteristic relaxation time of the alveolar. The timescale may describe a physical relaxation process greater effective than the values of and separately, such as unbinding rates or binding of molecular structures. The model parameterization (, , , , and ) was estimated by fitting to plotted stress–strain data inside the alveolar by experimental measurements. These experiments were produced from uniaxial tension measurements of the alveolar.

Discussion

Figure 3 in the graphical abstract shows step-by-step modeling of the respiratory system, including point cloud, CFD, FSI, particle deposition by discrete phase model (DPM) method, and final validation of the results using the particle image velocimetry (PIV) and stereolithography (SLA) method (Mortazavy Beni et al. 2019a, b, 2021; Mortazavi et al. 2020).

Fig. 3

The graphical abstract shows step-by-step modeling of the respiratory system according to previous literature

Respiratory diseases are associated with changes in respiratory viscoelastic behavior. The variation in ECM composition results in an alteration in the strength of the tissue. This variation was simulated by changing the volume fraction of ECM constituents in the FSI and MDM models. Changing the volume fraction of the pulmonary extracellular matrix (ECM) constituents leads to changes in elastic modulus (storage modulus) and the viscous modulus (loss modulus) of lung tissue. Therefore, in the biomedical computational methods where the model wall is considered flexible, the viscoelastic properties of the texture must be considered correctly. In general, in Fig. 4, the relationship between ECM properties changes and alterations in the mechanical characteristics of the lung tissue strip is well understood.

Fig. 4

The relationship between ECM properties changes and alterations in the mechanical characteristics of the lung tissue strip. a The extracellular matrix in health condition. b Spring-dashpot network model (infinite SLS) in health condition. c The extracellular matrix in disease. d Spring-dashpot network model (finite SLS) in disease

In emphysema, the divisions between the alveoli break down, producing larger lung spaces. This destruction of lung tissue reduces the springiness of the lungs. The lungs become more compliant and a small change in pressure produces a larger than normal change in the volume. While at first glance this would appear to make it easier to breathe, the opposite is true. Much of the work of breathing is done in overcoming the resistance of the airways. In emphysema, the airway resistance increases greatly. Also, increasing our understanding of these mechanical properties is essential for understanding lung function and feedback between FSI-MDM and lung tissue. This viscoelastic network can account for variations in the major ECM constituents, including elastin, collagen, PGs, and the surrounding matrix. The model demonstrates individual alveolar mechanical behaviors and links this into a network like the generalized Kelvin or Maxwell model. Relaxation and creep modulus simulations were performed to describe the bulk tissue strip mechanical manners. Reduction in the volume fraction of elastin, collagen, and PGs resulted in a decrease in the stiffness of the tissue. Also, in an experimental investigation (Yuan et al. 2000), it was illustrated that after digesting collagen, the stiffness of the tissue decreased. However, the interaction between collagen and elastin fibers is mediated via the surrounding matrix. Also, proof of a possible cross-link between these collagen and elastin fibers has been reported previously (Muiznieks and Keeley 2013). Hence, removing collagen will change the elastin structure and consequent response leading to a reduction in the stiffness of the entire tissue. On the other side, changes in body tissue properties can lead to changes in the behavioral patterns of various organs (Alaodolehei et al. 2020). However, the ECM state can be recognized by studying viscosity behavior; Hence, a more complex understanding of the respiratory system requires a profound insight into the fluidity and solidity of cells via the characterization of the ECM viscoelasticity. A viscoelastic model of the ECM was proposed with discrete differential equations after obtaining information about lung alveolar properties.

Conclusion

Increasing our understanding of mechanical properties is essential for understanding airway wall function and feedback between FSI-MDM and lung tissue. In this regard, considering a viscoelastic airway wall during transient airflow responding to the pressure pulse from the lung has a special attitude for future perspectives. The ECM contains a large part of the lung parenchymal cells. Its included three essential components: elastin and collagen within a proteoglycan (PG) viscoelastic network. Elastin provides the lung's elasticity property, a necessity for normal breathing, collagen prepares structural support and strength, and PGs give stability and cushioning within tissue loading. Respiratory diseases are associated with changes in the ECM composition. The variation in ECM composition results in an alteration in the strength of the tissue. For example, in emphysema the divisions between the alveoli break down, producing larger lung spaces. This destruction of lung tissue reduces the springiness of the lungs.