key: cord-0816963-bhp2d6z5
authors: Rajaraman, Prathish K.; Choi, Jiwoong; Hoffman, Eric A.; O'Shaughnessy, Patrick T.; Choi, Sanghun; Delvadia, Renishkumar; Babiskin, Andrew; Walenga, Ross; Lin, Ching-Long
title: Transport and deposition of hygroscopic particles in asthmatic subjects with and without airway narrowing
date: 2020-04-28
journal: J Aerosol Sci
DOI: 10.1016/j.jaerosci.2020.105581
sha: 0dface3b990e722daf092f9331a3b4814192dfe7
doc_id: 816963
cord_uid: bhp2d6z5

Abstract This study numerically investigates the effect of hygroscopicity on transport and deposition of particles in severe asthmatic lungs with distinct airway structures. The study human subjects were selected from two imaging-based severe asthmatic clusters with one characterized by non-constricted airways and the other by constricted airways in the lower left lobe (LLL). We compared the deposition fractions of sodium chloride (NaCl) particles with a range of aerodynamic diameters (1–8 μm) in cluster archetypes under conditions with and without hygroscopic growth. The temperature and water vapor distributions in the airways were simulated with an airway wall boundary condition that accounts for variable temperature and water vapor evaporation at the interface between the lumen and the airway surface liquid layer. On average, the deposition fraction increased by about 6% due to hygroscopic particle growth in the cluster subjects with constricted airways, while it increased by only about 0.5% in those with non-constricted airways. The effect of particle growth was most significant for particles with an initial diameter of 2 μm in the cluster subjects with constricted airways. The effect diminished with increasing particle size, especially for particles with an initial diameter larger than 4 μm. This suggests the necessity to differentiate asthmatic subjects by cluster in engineering the aerosol size for tailored treatment. Specifically, the treatment of severe asthmatic subjects who have constricted airways with inhalation aerosols may need submicron-sized hygroscopic particles to compensate for particle growth, if one targets for delivering to the peripheral region. These results could potentially inform the choice of particle size for inhalational drug delivery in a cluster-specific manner.

Inhalable aerosolized pharmaceutical drugs are commonly prescribed by medical doctors to be delivered throughout the respiratory tract for treating asthmatic patients. These drug aerosols experience mass change due to the interaction of aerosol and water vapor in the respiratory tract, which is a complex process because it depends upon local flow and particle characteristics, thermodynamics, airway structure, and lung function (Finlay, 2001; Longest & Kleinstreuer, 2005; Martonen, 1982) . A better understanding of the aerosol-water vapor interaction in the airways could potentially help the design of effective inhalers that engineer the aerosol size for targeted drug delivery (Byron, 2004; Longest, Tian, Walenga, & Hindle, 2012; Walenga & Longest, 2016) .

Recent studies have focused on the shrinkage or growth of particles for various airway conditions and its impact on the deposition of drugs in the airways (Chen, Feng, Zhong, & Kleinstreuer, 2017; Chen, Kleinstreuer, Zhong, Feng, & Zhou, 2018; Longest, McLeskey, & Hindle, 2010) . However, tailoring efficient drug delivery methods and devices for asthmatic patients poses numerous experimental challenges due to inaccuracy in measuring and quantifying the regional deposition of particles (dosage) in the airways (Borgström, Bengtsson, Derom, & Pauwels, 2000) . Although several computational fluid dynamics (CFD) studies have been performed for particle growth in the mouth-throat and airway model, very few have focused on the effect of subject-specific airways with realistic boundary conditions (Chen et al., 2017; Chen, Ma, Zhong, Sun, & Zhou, 2019; Winkler-Heil, Ferron, & Hofmann, 2014; . The change in the size of pharmaceutical particles is due to evaporation and condensation (hygroscopicity) when traveling in the respiratory tract (Broday & Georgopoulos, 2001; Finlay, 2001) . During the hygroscopic process of a sodium chloride (NaCl) particle, water vapor would be absorbed from surrounding air to form a film of saturated NaCl solution that covers a NaCl crystal (Asgharian, 2004; Broday & Georgopoulos, 2001; Chen et al., 2017) . As the hygroscopic process continues, the mass of solid NaCl decreases and the saturated NaCl solution increases until solid NaCl is fully dissolved (Asgharian, 2004; Broday & Georgopoulos, 2001; Chen et al., 2017) . Several previous works on the growth of hygroscopic particles in the human airways showed that an increase in particle size enhances the particle deposition throughout the airways (Chen et al., 2019; Ferron, Kreyling, & Haider, 1988; Ferron, Oberdörster, & Henneberg, 1989; Hofmann, 2011; Martonen, 1982) . For example, Ferron et al. (1989) found that for particles with an initial diameter of 0.5-2 μm in a constant relative humidity (RH) environment, the hygroscopic process could significantly alter the particle size by a factor of 2 to 4. Chen et al. (2017) tested the hygroscopic growth of multi-component particles, comprising water, ethanol, NaCl and fluorescein, in a generalized mouth-throat model with various boundary conditions at the inlet and boundary walls, such as volumetric flow rate and wall wetness. They found that for particles with an initial diameter of 2.5-7.5 μm a wet wall (RH = 99.5%) condition increases the deposition efficiency by 2.1%, as compared with 0.5% obtained with a dry airways wall condition (RH = 0%). In a numerical study of hygroscopic particle growth of NaCl, Chen et al. (2019) reported that a constant RH of 99.5% at the airway wall overpredicted the final particle diameter by 18% in the airways of generation 3-6, indicating the importance of this boundary condition.

Cluster analysis allows for grouping subjects with similar imaging-based structural, functional metrics and clinical characteristics, such as pulmonary test functions (PFTs), body mass index (BMI), activity limitation and exacerbations, facilitating bridging individual and population-based studies. Thus, the notion of subject-specificity for tailoring drug delivery could potentially be realized through cluster-specificity. For example, Choi et al. (2017) identified four asthmatic clusters using computed tomography (CT) imaging-based variables and established their associations with clinical characteristics. Later, J. Choi et al. (2019) applied a CFD model to examine air flow and deposition of non-growth particles within archetypes representative of the four imaging-based clusters. Although subjects in cluster 3 (without airway constriction) and cluster 4 (with airway constriction) were characteristic of severe asthma, they exhibited different particle deposition patterns with increasing initial particle size. Specifically, cluster 4 subjects had an increase in the deposition fraction in the left lower lobe (LLL) for particles with a diameter larger than 2 μm.

This study focuses on the effect of particle growth on particle transport and deposition in severe asthmatic subjects with and without airway constriction, which is an extension of previous work (J. Choi et al., 2019) . A CFD model was implemented to predict the transport and deposition of dry NaCl particles in archetypes selected from two severe asthmatic clusters.

Cluster 3 (C3) subjects had non-constricted airways, while cluster 4 (C4) subjects had constricted airways. In this study, only the growth of NaCl particles was considered due to the high affinity of NaCl particles to absorb moisture from the surrounding air compared to pharmaceutical drugs.

The effects of particle growth for a range of dry initial particle sizes in each cluster were investigated on the deposition pattern and final particle sizes in the airways, specifically in LLL, where airway narrowing exists.

Computed tomography (CT) lung images of six severe asthmatic subjects at total lung capacity (TLC) and functional residual capacity (FRC) were acquired at the University of Pittsburgh and the University of Wisconsin with approval from respective institutional review boards, as part of the Severe Asthma Research Program (SARP). These subjects are archetypes representative of the two severe asthma clusters. Females with non-constricted airways and reversible lung function dominated C3. In contrast, C4 was dominated by older males with constricted airways and persistent reduced lung function (Choi et al., 2017) . For each cluster, only three subjects, denoted by (i), (ii) and (iii) in Table I, were selected for CFD analysis owing to the high computational cost for simulation of turbulent flows. Their features are summarized in Table I , such as sex, age, reduced * ( ), and reduced ( ) as marked by arrows ↓ for C4. * ( ) is the average of normalized airway hydraulic diameter in the LLL. is the Jacobian determinant defined as the ratio of a local lung volume at TLC over that of the same tissue mass at FRC and ( ) is the average Jacobian determinant over the whole lung.

We employed a large eddy simulation (LES) technique that solves the following filtered governing equations (Lin, Tawhai, McLennan, & Hoffman, 2007; Wu, Tawhai, Hoffman, & Lin, 2014) :

-

where , andare the filtered velocity vector, temperature and water vapor concentration, respectively. is the fluid density, ! is the filtered pressure, % is the fluid viscosity, and % & is the subgrid-scale (SGS) eddy viscosity proposed by Vreman (2004 (8) 7 ) of water vapor in the air were fixed at 0.9 . The buoyancy force is modeled by the Boussinesq approximation as " = −9 ′ . ; , where 9 is the coefficient for volumetric expansion of air, ′ . is the modified temperature of dry air accounting for moist air and ; = (0,0, −;), where ; is the gravitational acceleration. The Reynolds number in the trachea (Re trachea ), based on the average tracheal diameter and the maximum volumetric flow rate at peak inspiration = 0> ( /@AB), varied from 2,742 to 4,881 ( Table Ⅰ) . The = 0> was derived from subject-specific inspiratory capacity (J. Choi et al., 2019).

The three-dimensional (3D) segmented human airway geometry and one-dimensional (1D) skeleton were obtained from CT images at TLC for each cluster subject using Apollo software (VIDA Diagnostics). CT scans are acquired when the study subject lies on the imaging bed in a supine position. The gravitational direction for the CFD analysis is taken as the axis of rotation of the CT scanner, parallel to the imaging bed. Thus, the gravitational direction is not exactly parallel to the trachea. The 3D segmented airway geometry and skeleton were then input into an in-house code to automatically construct a subject-specific 3D-1D coupled airway geometrical model, which contains anatomical information and meshes for CFD (Miyawaki, Tawhai, Hoffman, Wenzel, & Lin, 2017) . The CT images, however, did not include the upper airways. To mimic the effects of the upper airways, we employed a laryngeal model to create a glottal constriction in the trachea that could generate a turbulent jet on inspiration .

The meshes generated had approximately 8 to 14 million tetrahedral elements, as shown in Figure 1 . The meshes used in this study was generated using Gmsh (Geuzaine & Remacle, 2009 ) and Tetgen (Si, 2015) , which used the 3D subject specific airway geometry with the 1D airway tree information. The 1D airway tree for all the subject was generated using a volume filling technique (Tawhai, Pullan, & Hunter, 2000) . The mean branching angle is fairly constant (37 ο ) and the branch length varied from 2 mm to 28 mm for Horsfield order 1-20 comprising sublobar bronchi and bronchioles. The 3D mesh grid size was computed using the branch diameters and flow rates in order to preserve the grid size in terms of wall units in every branch . The grid size in wall units D E was computed based on the following equation.

Where ∆G H is the branch grid size, 6 J,H is the flow rate ratio of branch to trachea and H is the branch diameter. A non-uniform mesh was generated in each branch to capture the turbulent effect accurately, for more detail on the meshing algorithm implemented see .

In order to obtain a subject-specific flow distribution, we employed a mass-preserving image registration method that matches CT scans at TLC and FRC to estimate the air volume change at each terminal bronchiole of the 1D tree (Yin, Choi, Hoffman, Tawhai, & Lin, 2010) .

The air volume change data were used to derive local compliance coefficients for a 1D resistance model. Then, dynamic flow rate distributions at terminal bronchioles were derived from a 1D resistance simulation, which were then used to generate CFD flow rate boundary conditions for CT resolved airways (S. Choi et al., 2019) . Contrary to previous numerical studies that used a fixed volumetric flow rate (= 0> ) at peak inhalation (Chen et al., 2017; Chen et al., 2019; Feng, Kleinstreuer, Castro, & Rostami, 2016; Zhang, Kleinstreuer, Kim, & Cheng, 2004) , we employed the subject-specific maximum volumetric flow rate (= 0> ) derived from the CT images, with inspiratory capacity (IC) shown in Table Ⅰ . We fixed the peak inhalation time at 1.25 s to mimic the characteristic inhalation waveform for slow and deep breathing using an MDI with the inhalation period of 5.0 s for all subjects. The waveform was proposed by Longest et al. (2012) based on package insert instructions. The = 0> value is essential when studying hygroscopic particle growth since the flow rate determines the amount of time for particle-water vapor interaction. At the airway wall, a no-slip boundary condition was imposed.

A unique feature of our approach lies in the modeling of temperature and water vapor concentration in the airways. Most of the existing studies assumed a wet wall condition where the temperature (37Ⅰ C) and relative humidity (RH = 99.5%) were fixed at the airways wall (Asgharian, 2004; Broday & Georgopoulos, 2001; Ferron et al., 1988; Longest et al., 2012) . A fixed airway wall temperature or using simple function of distance from mouth may lead to inaccurate prediction of temperature distribution within the airway since during inspiration the airway wall is cooled due to evaporation and heat transfer (Winkler-Heil et al., 2014; . The wall boundary conditions for Equations (3) and (4) and humidity are transferred to the air from the airway wall until isothermal saturation is reached (Wu et al., 2014) . For Equations (3) and (4), a Dirichlet boundary condition of 303 K and RH = 35% was imposed at the entrance of the trachea. At the end of the terminal bronchioles, a Neumann condition was imposed (Wu et al., 2014) . Once the velocity, temperature, and water vapor concentration (absolute humidity) fields were obtained, the particle transport and growth models described in the following sections were solved.

The trajectories for spherical particles were computed using a Lagrangian particletracking algorithm. The particle velocities computed here were based on the assumption that particles do not affect airflow, gas phase temperature, and water vapor concentration, since a dilute regime is considered here:

where * is the particle velocity, " M is the drag force per unit mass, and * is the particle density (Lambert, O'Shaughnessy, Tawhai, Hoffman, & Lin, 2011) . The total number of particles considered for this study was 200,000. Both Brownian and Saffman forces, which are only important for submicron particles (Li & Ahmadi, 1992) , were omitted due to the particle diameters considered in this study (1-8μm). They were released at the beginning of the simulation uniformly in a cross-section at the entrance above the glottal constriction. The condition for particle deposition to occur is when the distance between the particle center and the airway wall is less than the particle radius. Equation (6) was solved for a single inspiration (1.25 seconds) to mimic slow and deep inhalation of a metered dose inhaler.

Inhaled aerosol particles experience mass/size change in the human respiratory tract, and the mass change occurs due to evaporation/condensation of each component (N) in a particle until equilibrium is reached (Chen et al., 2017; Chen et al., 2018; Longest & Kleinstreuer, 2005; Worth Longest & Hindle, 2011; Worth Longest & Xi, 2008) . In this study, we considered the component of water vapor only. The rate of particle mass change when traveling in the airways is as follows:

where @ * is the mass of the particle, B Q is the average mass flux of evaporable/condensable, of water vapor on the particle surface ( ) due to the presence of component (N) and k is the total number of components (Chen et al., 2017) . More often, the particle mass change in Equation (7) can be reformulated in terms of particle diameter $L * ' given by (Longest & Kleinstreuer, 2005) :

In Equation (8), the average mass flux is defined as:

where Z Q is the diffusion coefficient of component N in air and 8ℎ is the Sherwood number calculated using Equation (10), which is valid for Re c < 400 (Clift, Grace, & Weber, 2005) .

where 8) is the Schmidt Number and fN * is the particle Reynolds number. -[ in Equation (9) is the Fuchs-Knudsen number correction, defined below where i is the gas mean free path, and j [ is the mass accommodation coefficient (Chen et al., 2017) :

In Equation (9), ] Q,S_1 and ] Q,^ are the mass fractions of each component (N) on the surface of the particle and in the gas phase far from the particle, respectively. The mass fraction in the gas phase $] Q,^' is determined by Equation (4), which is computed by taking the average of the nodal value and interpolating to the location, where the particle currently resides. The mass fraction of water vapor on the particle surface can be computed from the equation below:

where q Q is the activity coefficient to account for the deviation from Raoult's Law, G Q is the mole fraction, f Q is the universal gas constant, * is the particle temperature, 5 .Q,Sr7 $ * ' is the saturation pressure and + Q is the Kelvin effect, which is the correction factor included here to account for curved surfaces (Finlay, 2001; Longest & Kleinstreuer, 2005) .

In Equation (13), the heat transfer of a particle to surrounding air is derived from the energy balance equation when the heat is gained or lost by the particle due to phase change and heat exchange between * and the surrounding air temperature ( r ) in the respiratory tract (Asgharian, Rostami, Price, & Pithawalla, 2018; Chen et al., 2017; Finlay, 2001) . In Equation (13), i t u is the thermal conductivity of air-vapor mixture. ) w is the particle heat capacity, and Q is the heat of vaporization/condensation of each particle component (Chen et al., 2018) . The

Nusselt number (v ) is expressed below, where 56 is the Prandtl number (Chen et al., 2017) .

Quantification of the number of particles depositing on the respiratory tract wall requires the implementation of appropriate conditions at the airways wall. The airway wall here is approximated as a sink, which assumes nearly instantaneous reaction kinetics between the wall and particles that the particle growth is stopped instantly (Finlay, 2001) . The particle deposition in the airways can be quantified in terms of deposition fraction (DF), as given in Equation (15).

The total number of particles depositing in the whole lung or lobe $v *,wQ* ' is divided by the total number of particles entering $ v *,Qz7 ' at the top of the laryngeal model (supraglottal region). The deposition density (DD), as expressed in Equation (16), where R S is the unit surface area of the airways. The effect of particle growth to non-growth particle was accessed by computing the relative difference in DF, as shown in Equation (17). Throughout this study, the average DF between cluster subjects was computed. Because there are three subjects in each cluster, the average values are presented for analysis.

For validation of the above particle growth model, we compared our model results with the experimental data from a previous study (W. Li, Montassier, & Hopke, 1992) . In W. , the size of the NaCl particles was sorted using a differential mobility analyzer (DMA).

The particles were then allowed to flow into a wetted-wall tube reactor where the hygroscopic growth occurred with the air stream fixed at RH > 99%. Once the particles reached an equilibrium state, the aerosol diameter was measured with a second DMA associated with a condensation nucleus counter. For any given initial particle diameter, the experiment was repeated multiple times to obtain a mean diameter of the hygroscopically grown particle. In Table Ⅰ , the parameters used were obtained from the study (Chen et al., 2017) . We assumed a steady-state condition by fixing the temperature and humidity of the air stream entering the simulated tube reactor to 37 o C, RH = 99% and 22 o C, RH = 99.5% (W. . The numerical simulation used dry NaCl particles with initial diameters ranging from 10 to 200 nm.

The model used here assumed a rapid mixing of dissolved salt in the droplet so that a concentration gradient within the particle is ignored, which is a reasonable assumption for NaCl when modeling hygroscopic particle growth (Asgharian, 2004; Broday & Georgopoulos, 2001; Chen et al., 2018) . In Equation (12) (Chen et al., 2017; Finlay, 2001) .

In Figure 2 , we compared the predicted particle growth ratios with the experimental data (W. for dry NaCl particles with initial diameters ranging from 10 to 200 nm at two different conditions. The growth ratio here was calculated by taking the equilibrium particle diameter $L * ' and dividing by the initial particle diameter $L w1" '. The model prediction at both conditions shows an increase of growth ratio with increasing initial particle diameter. The predicted growth ratios here are for a fixed RH = 99% (solid black circle) and RH = 99.5% (solid black square). Previous studies showed that a small change in relative humidity (RH > 98%) could significantly impact the growth ratio and the time to reach the equilibrium size (Broday & Georgopoulos, 2001; Haddrell et al., 2017; Winkler-Heil et al., 2014) . The discrepancies between the model prediction and experimental data are due to several factors. First, the growth ratio shown in Figure 2 heavily depends on the RH, which was measured by W. with an accuracy of 1% as RH approaches 100%. As shown in Figure 2 , a change in RH by 0.5% increased the growth ratio by approximately ~1.0, when the initial particle diameter was fixed at 200 nm. Second, the sensitivity of particle growth model discussed above to thermophysical properties such as diffusion coefficient $ Z Q ', thermal conductivity $ i Z u ', and saturation vapor pressure $5 .Q,Sr7 ' are well-studied, indicating a difference in growth ratio about ~1.45 using thermophysical properties from various studies (Broday & Georgopoulos, 2001; Chen et al., 2017) . Although there exists a mismatch between the model prediction and experimental data, the agreement between them could be further improved by adjusting the RH in the growth model by ± 0.2-0.5%, depending on initial particle diameter.

In Figure 3 , we compared the iso-surfaces of airflow speed for subject C3(Ⅰ) in ( and (f), which is about 0.039 kg/m 3 for the C3 subject and 0.035 kg/m 3 for the C4 subject. A slight decrease in local water vapor concentration is associated with an increase in airflow speed and a decrease in temperature. The water vapor concentration in the airways has a greater impact on the particle trajectory since different water vapor concentration may lead to different particle size, which subsequently affects the deposition site in the airways due to size change and particle momentum. Figure 4 compares the average DF values for the C3 and C4 subjects with and without (constant particle diameter) hygroscopic particle growth for various initial particle diameters for a simulation time of 1.25 s. The average DF increases as the initial particle size varied from 1 to 8 μm since the particles that would be advected deeper into the airways now deposited on the airway wall. This is true for both cases with and without hygroscopic particle growth. As expected, with particle growth, the average DF increases. The difference in DFs between with and without growth is most evident for particles with an initial diameter of 4 μm, viz. about 3% in the C3 subjects and 6% in the C4 subjects. Comparing the average DF between the C3 subjects (Figure 4(a) ) and the C4 subjects (Figure 4 (b)) for all initial particle diameters, on average, the C4 subjects have higher DF values because of constriction in the C4 subjects. In addition, the average growth ratio decreases from about 2.0 to 1.15, with increasing initial particle size from 1 to 8 μm, which was computed for all particles that were deposited or advected to the outlet. In Figure 4 , the decrease in growth ratio was observed due to a lower surface to volume ratio, which slowed the response of particles in a humid environment. In order to test the effect of particle release time, we compared the deposition fractions for subjects with constriction. We released monodispersed particles with initial diameters ranging from 1 to 8 μm at 0.3 s. Our findings showed an increase in deposition fraction about 3.5%-5.0% compared to releasing all the particles at the initial time. The increase in deposition fraction is due to the turbulent intensity experienced by particles at a later time. Our observation suggests that releasing a fixed number of particles at each time, will likely increase the deposition for particles released closer to the peak inhalation time. The increase in deposition occurs since higher velocity leads to more inertia.

There are several factors influencing the particle growth in the airways. First, the volumetric flow rate for each subject in this study is subject-specific as described above, which influences the particle-water vapor interaction time, since a larger = 0> leads to a smaller growth ratio for a micron-sized particle. Second, the sensitivity of the particle growth model to physical and chemical properties of water was well studied by (Chen et al., 2017) . For example, Chen et al. (2017) showed that for a small change of 50 Pa in saturated water vapor pressure, the estimated particle growth ratio was about ~1.45 fold larger. As shown in Figure ( 3) ((c)-(f)) the variation in both water vapor concentration and temperature due to evaporation would have a greater influence on the growth curve ( Figure 4 ) compared to predictions modeled with constant temperature and water vapor concentration, as described in Equation (9) where ] Q,^ is computed.

A more recent study by Chen et al. (2019) concluded that the impact of variable RH in generalized airways of generation 3-6 can impact the final particle diameter up to 18%. This highlights the importance of accurate prediction of water vapor concentration. Furthermore, the influence of lung geometry also contributes to the growth ratio, since the variations in local airway diameter throughout the airway tree, as well as local velocity, lead to a distribution of growth ratio in each generation.

The difference between the predicted average DF for both C3 and C4 subjects shows the importance of including hygroscopic particle growth model, accurate wall temperature and water vapor concentration, which affect the local distribution of particles/drugs within the airways (dosage). Many studies have shown that when modeling particles in the airways for a single breath simulation, submicron particles achieve equilibrium size. These studies are either experiments conducted in a wet wall tube or computational models that used a generalized airway geometry with constant relative humidity and temperature throughout the airways. A comparison of the average DF values shown in Figure 4 with the study by Asgharian (2004) highlights several differences between them. In Asgharian (2004), the hygroscopic particle growth of NaCl was compared for a breathing frequency of 15 breaths per minute with a typical oral breathing tidal volume of 625 ml. The DF in the tracheobronchial region for an initial particle diameter of 2 μm was approximately 36% with a symmetric lung model compared to the current model with an average DF of 13% in the C3 subjects and 23% in the C4 subjects. The average DF value presented in this study is considerably lower than the value reported by Asgharian (2004) .

There exist several key differences between the model presented by Asgharian (2004) and the current model. First, the current model accounts for the movement of the particle in the fluid by correcting the mass flux term in Equation (9) with 8ℎ , whereas Asgharian (2004) excluded 8ℎ. Second, their model assumed constant temperature and RH in the airways, which allowed particles to reach equilibrium in a short time contrary to the current model. Third, the model implemented a symmetric airway geometry, while the current airway geometries are subject-specific reconstructed from CT images. In the study by Hindle and Longest (2010) , the authors performed a comparison of experimental and numerical growths of albuterol sulfate particles in an idealized airway geometry for initial particle diameters of 0.56 μm and 0.9 μm.

They predicted the DF of 8.2% for 0.9 μm initial particles in a geometric lung model including the mouth-throat region. In comparison, the DF shown in Figure 4 for an initial particle diameter of 1 μm is approximately 10% in the C3 subjects and 16% in the C4 subjects. They were in a good agreement for the C3 subjects in spite of the difference in particles, albuterol sulfate particles in Hindle and Longest (2010) vs. NaCl particles in this study.

To examine the impact of the hygroscopic particle growth on DF in local regions, we computed the relative difference in average DF ("enhancement" due to hygroscopic growth)

using Equation (17) in the entire upper (RUL, LUL) and lower lobes (RML, RLL, LLL) for a range of initial particle diameters, as displayed in Figure 5 for (a) upper and (b) lower lobes. That is, the relative difference in average DF between the cases with and without growth can be interpreted as the deposition enhancement due to hygroscopic particle growth when compared to the case without hygroscopic particle growth. In the upper lobes ( Figure 5 (a) ), for initial 1 μm particles in the C3 subjects, the relative difference between with and without particle growth shows a small ∆DF (~5%). The 1 μm particle in the C4 subjects (square) shows an increase about ~30% relative to non-growth particles as a result of water condensation while in the upper lobes. The most substantial difference between growth and non-growth is for particles with an initial size of 2 μm. For particles with initial diameter larger than 2 μm, ∆DF decreases. In the lower lobes, the difference is the highest for the C4 subjects ( Figure 5 (b)) with an initial particle diameter of 2 μm (~70%), indicating that hygroscopic particle growth is more sensitive in lobes characterized by airway narrowing.

Excluding the growth model for 2 μm particles, the mechanism that dominates for particle deposition was mostly due to sedimentation. With the growth model, 2 μm particles grow due to accumulation of water and the dominating mechanism for deposition becomes impaction. In comparison with Figure 4 , the difference in average DF for initial 4 μm particles appears to be larger than those of 2 μm. This is attributed to higher deposition in the carina region, which was excluded when computing the ∆DF in lobar regions ( Figure 5 ). In summary, the growth model has more impact on the C4 subjects, especially in the lower lobes where airway narrowing exists. Figure 6 shows the final particle diameter and deposition location after hygroscopic growth for representative cluster subjects C3(ii) and C4(iii) in Table I at peak inhalation time of 1.25 seconds. The initial particle diameter was fixed at 2 μm and the particle growth predicted a higher final particle diameter in the trachea region (more red particles), varying from 3.13 to 3.50

μm. This is due to the turbulent inspiratory flow that increases heat and mass transfer. Once the flow passes the carina region, the flow is split into two streams at the first bifurcation, where particles tend to travel in a circumferential manner to deposit at the bottom of the first bifurcation.

Particles that are transported deeper into the lungs till the next bifurcation where more particles are deposited at the ridges (via impaction) and active particles are further split into two groups, which travel into the upper and lower airways. In the lower airways (see the boxed area in Figure   6 ), DF increases with increasing initial particle diameter due to impaction. In Figure 6 (a), the boxed region focuses on the LLL region where the airways are classified as having a nonconstricted airway diameter using an imaging-based clustering algorithm. The boxed region shows that the particles deposited with final particle diameter ranging from 2.88 to 3.25 μm have more time for particle-water vapor interaction. In Figure 6 (b), we can see an increase in particle deposition in the boxed region due to airway narrowing with the final particle diameter ranging from 2.75 to 3.13 μm.

In Figure 7 , the DD is shown for two representative cluster subjects. As expected, DD is consistently higher in the C4 subjects compared to the C3 subjects at the constriction. For the C3 subject shown in Figure 7 (a), the DD in the boxed region appears in small areas, and mostly on the lower airway wall due to non-constricted airways and associated lower airflow speed in this region. As for the C4 subject (Figure 7 (b)) the effect of particle growth at the constriction is more evident because of an increase in particle mass/size and airflow speed and a decrease in surface area. The particles that escaped the constriction traveled deeper into the airways and continued to grow compared to the particles deposited in the boxed region.

Incorporating the effects of particle growth into numerical models provides valuable information for targeted drug delivery in severe asthmatic cluster subjects. The impact of particle growth shown in Figure 6 indicates that the constriction in the LLL for the C4 subjects has considerable effect in terms of particle deposition with an initial particle diameter of 2 μm. To reduce this effect in C4 subjects, the initial particle diameter required may be smaller, enabling the delivery of more particles to the alveolar region even accounting for particle growth.

Previous studies have suggested that particles with an initial diameter < 1.0 μm are sufficient to escape the upper airways and grow to 2.5 μm when entering the alveolar region (Worth Longest & Hindle, 2011) . Based on the current study, submicron particles would reduce the DF in the LLL with constricted airways, since the final particle size would be less than 1 μm even with growth and more particles could be advected to the alveolar region. In C3 subjects with nonconstricted airways, limiting the initial particle diameter to < 4.0 μm, the overall DF in the airways was less than 12%. The model presented here predicts the aerodynamic and thermodynamic mechanisms of NaCl particles. In the future, the current model can be extended for particles with multi-components representing drug aerosols.

The inter-subject variabilities of deposition fraction in human airways are very complex.

We aimed to shed light on inter-cluster analysis with inter-subject variabilities, where clustering characterized discriminant structural and functional features in subgroups. The variability in airway geometry (length, diameter), breathing pattern, and regional ventilation, also influences the temperature and humidity distribution (Wu et al., 2015) . Even though more subjects are preferable for CFD analysis to obtain a complete statistical analysis, we limited the total number of subjects to six due to computational cost. Specifically, we focused this study on C3 and C4 subjects, which were distinguished by distinct airway structures and DF (J. Choi et al., 2019) . In Table Ⅰ , the means of important imaging-based metrics, CFD variables, and deposition fractions in LLL for the growth model were compared between these cluster subjects using biserial correlation analysis. The analysis showed that the cluster memberships of selected C3 and C4 subjects were discriminated by * ( ) with a negative correlation of r = -0.79 (p = 0.058).

Additionally, the deposition fractions in LLL, {( ), were greater in C4 subjects than C3 subjects, especially for large initial particle sizes of 4μm (r = 0.80, p = 0.056) and 8μm (r = 0.82, p = 0.048). Overall, the functional and structural variables (J and * ) were smaller in C4 subjects than C3 subjects and inversely correlated with {( ). The result is consistent with previous findings without the particle growth model (J. Choi et al., 2019) . {( ) had a stronger correlation with * ( ) than ( ) and = 0> . * ( ) had significant negative correlations that became stronger with increasing particle size. ( ) showed significant correlations only at large particles of 4μm (r = -0.44, p = 0.045) and 8μm (r = -0.47, p = 0.004). The correlations with the peak inhalation trachea flow rate ( = 0> ) were significant at all sizes, but became weaker with increasing particle size.

The presence of airway narrowing increases the inter-cluster variability of water loss at the airway wall, as shown in Figure 3 (e) and (f). Thus, this highlights the importance of using subject-specific geometry and ventilation. In a previous study, Wu et al. (2014) proposed dimensionless correlations based on local flow rate (fN) , local airways diameter ( r ) and trachea diameter ( 7 ). n increases. The altered local mass and heat transfer in the fluid phase are linked to the particle phase via equations (9) and (13). The effect of the constriction can be seen in Figure 6 . For the C4 subjects, the local flow velocities increases in the constriction (Figure 3) , the particles entering the constriction experienced more inertia, which causes the particles to deposit at the constriction. Even though higher Nu LLL and Sh LLL are observed at the constriction (C4), the particle-water vapor interaction time is lower in the constriction, which indicates smaller deposited particle diameter (Figure 6(b) ). This effect does not occur in the C3 subjects due to the lower flow velocities at LLL. For C3 subjects, the particles entering LLL gains mass via condensation, and the final deposited particle diameters are larger, see Figure 6 (a). Thus, this highlights the importance of predicting water evaporation at the airways wall, since the mass of water transferred to the particle surface depends upon the humidity of the surrounding air. In Table Ⅰ , Pearson's correlation analysis for Nu LLL and Sh LLL with deposition fractions in LLL are shown, which indicates a strong positive correlation for deposition fraction for 1 μm case (r = 0.98, p < 0.001 ) and 2 μm (r = 0.96, p < 0.001). The initial particle diameter increases the correlation (Nu LLL and Sh LLL ) decreases (Table Ⅰ) , which is due to lower growth rates for larger particle diameters, as shown in Figure 4 . The trend is observed due to the simulation time of 1.25 s (single inhalation), where larger particles require a longer time to reach equilibrium. Airway constriction found in the C4 subjects leads to lower particle delivery to the peripheral airways, which may limit the performance of asthma inhalers.

The present work differs from other studies in three aspects. First, the subject-specific airway geometry was reconstructed from CT images, and the associated subject-specific inhalation capacity was computed using image registration. Furthermore, the boundary condition imposed at the terminal CT-resolved 3D airways for the 3D CFD simulation was computed based on a 1D resistance model, which could reflect the effect of local airway constriction in asthmatic subjects. Second, we employed a two-layer wall boundary condition to account for evaporation at the ASL-lumen interface. Third, we focused on investigating the hygroscopic effect in two distinct severe asthmatic subgroups, one without airway constriction and the other with constriction. It is noted that the cluster-specific structural and functional features are not detectable by clinical spirometry. These results could potentially inform the choice of particle size for inhalational drug delivery in a cluster-specific manner.

This study has four major limitations. First, the number of subjects in each cluster was limited to three due to high computational time (6 to 7 days per subject). Thus, a complete statistical analysis was not presented. The primary goal of this study was to examine the impact of particle growth in severe asthmatic cluster subjects with and without constriction in the LLL.

Second, we did not consider the effect of the mouth cavity due to the lack of CT image data . Third, a particularly interesting problem, which is not addressed in this study, is the effect of particle coagulation on DF (Asgharian et al., 2018) . Finally, the initial state of the particle was assumed to be dry at the inlet. Depending on the initial particle diameter and the composition of the particle, wet particles at the inlet may shrink or grow depending on the local humidity in the airways, and the change in particle diameter will impact the particle trajectory and deposition fraction in the airways.

The objective of the present work was to evaluate the effect of hygroscopic NaCl particle deposition fraction in CT-based subject-specific airway models from imaging-based severe asthmatic clusters. The impact of hygroscopic particle growth was compared with non-growth particles by computing the total deposition fraction values in each subject. The effect of particle growth on the DF in severe asthmatic subjects with non-constricted airway diameter was overall less than in subjects with constricted airways. The effect of hygroscopic particle growth model was evaluated by taking the differences in average DF for the cases with and without growth. For C4 subjects, the hygroscopic particle growth with a fixed initial particle diameter of 2 μm contributed more to the DF, which is approximately 3.5% higher in the LLL. The higher DF in the C4 subjects is due to the local conditions such as reduced airway diameter, velocity, temperature, and available water vapor. It is important to note that the choice of using either the growth model or no-growth model is much more important for C4 subjects than C3 cluster subjects. Predicting accurate temperature and water vapor concentration in the airways is crucial for computing hygroscopic particle growth and deposition fraction, especially in severe asthmatic subjects with airway narrowing. Utilizing cluster-specific imaging-based features to group severe asthmatic subjects into clusters such as airways constriction may help tailor clusterspecific hygroscopic drug treatment for improving inhalation drug delivery. The proposed CFD lung model with hygroscopic particle growth has potential applications for studies related to exposure to tobacco cigarettes, electronic cigarettes, environmental particulate matters, and virus-laden aerosolized ambient particles, such as COVID-19.

mention of trade names, commercial practices, or organization imply endorsement by the United States Government.

Eric A. Hoffman is a shareholder in VIDA diagnostics, a company that is commercializing lung image analysis software derived by the University of Iowa lung imaging group. He is also a member of the Siemens CT advisory board. 10, 20, 30, 50, 75, 100, 125, 150, 175, 200 Table . Summary of means (± SD) of significant imaging-based metrics, CFD variables and deposition fractions in LLL with particle growth for C3 and C4 subjects. The relationships between variables and clusters are computed using biserial correlation coefficient (r) with p < 0. Table . Pearson's correlation coefficients (p-value) between key variables and deposition fractions in LLL at four initial particle diameters. . Distributions of deposited particles color-coded by particle diameter for (a) C3(ii) and (b) C4 (iii) at 1.25 seconds for an initial particle diameter of 2 μm. The dashed boxed region highlights a segment of LLL (a) non-constriction and (b) constriction.

(a) (b) Figure 7 . The contours of the deposition density, defined as the number of deposited particles , per unit surface area (mm 2 ), for cluster subjects C3(ii) and C4(iii), in (a) and (b), with an initial particle diameter of 2 μm. The boxed region compares the number of particles deposition in the region (a) nonconstriction and (b) constriction.

A Model of Deposition of Hygroscopic Particles in the Human Lung

Regional deposition of inhaled aerosol constituents from Electronic Nicotine Delivery Systems (ENDS) in the respiratory tract

Variability in lung deposition of inhaled drug, within and between asthmatic patients, with a pMDI and a dry powder inhaler

Growth and Deposition of Hygroscopic Particulate Matter in the Human Lungs

Drug delivery devices: issues in drug development

Numerical investigation of the interaction, transport and deposition of multicomponent droplets in a simple mouth-throat model

Effects of thermal airflow and mucuslayer interaction on hygroscopic droplet deposition in a simple mouth-throat model

Numerical study of the effects of temperature and humidity on the transport and deposition of hygroscopic aerosols in a G3-G6 airway

Differences in Particle Deposition Between Members of Imaging-Based Asthma Clusters

Quantitative computed tomographic imaging-based clustering differentiates asthmatic subgroups with distinctive clinical phenotypes

1D network simulations for evaluating regional flow and pressure distributions in healthy and asthmatic human lungs

Bubbles, drops, and particles: Courier Corporation

Computational transport, phase change and deposition analysis of inhaled multicomponent droplet-vapor mixtures in an idealized human upper lung model

Inhalation of salt aerosol particles-II. growth and deposition in the human respiratory tract

Estimation of the deposition of aerosolized drugs in the human respiratory tract due to hygroscopic growth

The mechanics of inhaled pharmaceutical aerosols: an introduction

Gmsh: A 3-D finite element mesh generator with built-in pre-and post-processing facilities

Pulmonary aerosol delivery and the importance of growth dynamics

Evaluation of Enhanced Condensational Growth (ECG) for Controlled Respiratory Drug Delivery in a Mouth-Throat and Upper Tracheobronchial Model

Modelling inhaled particle deposition in the human lung-a review

Regional Deposition of Particles in an Image-Based Airway Model: Large-Eddy Simulation and Left-Right Lung Ventilation Asymmetry

Dispersion and deposition of spherical particles from point sources in a turbulent channel flow

A System to Measure the Hygroscopicity of Aerosol Particles

Characteristics of the turbulent laryngeal jet and its effect on airflow in the human intra-thoracic airways

Computational Models for Simulating Multicomponent Aerosol Evaporation in the Upper Respiratory Airways

Characterization of Nanoaerosol Size Change During Enhanced Condensational Growth

Comparing MDI and DPI Aerosol Deposition Using In Vitro Experiments and a New Stochastic Individual Path (SIP) Model of the Conducting Airways

Analytical model of hygroscopic particle behavior in human airways

Numerical simulations of aerosol delivery to the human lung with an idealized laryngeal model, image-based airway model, and automatic meshing algorithm

Automatic construction of subject-specific human airway geometry including trifurcations based on a CT-segmented airway skeleton and surface

TetGen, a Delaunay-based quality tetrahedral mesh generator

Generation of an Anatomically Based Three-Dimensional Model of the Conducting Airways

An eddy-viscosity subgrid-scale model for turbulent shear flow: Algebraic theory and applications

Current Inhalers Deliver Very Small Doses to the Lower Tracheobronchial Airways: Assessment of Healthy and Constricted Lungs

Calculation of hygroscopic particle deposition in the human lung

Numerical Model to Characterize the Size Increase of Combination Drug and Hygroscopic Excipient Nanoparticle Aerosols

Condensational Growth May Contribute to the Enhanced Deposition of Cigarette Smoke Particles in the Upper Respiratory Tract

A Numerical Study of Water Loss Rate Distributions in MDCT-Based Human Airway Models

A Numerical Study of Heat and Water Vapor Transfer in MDCT-Based Human Airway Models

Simulation of pulmonary air flow with a subject-specific boundary condition

Water Vapor Transport and Its Effects on the Deposition of Hygroscopic Droplets in a Human Upper Airway Model

Isotonic and Hypertonic Saline Droplet Deposition in a Human Upper Airway Model

Vaporizing Microdroplet Inhalation, Transport, and Deposition in a Human Upper Airway Model

This work was supported in part by FDA grant U01-FD005837, NIH grants U01-HL114494, and P30ES005605. We also thank the Extreme Science and Engineering Discovery Environment