Introduction

Inside of our nose is structurally and physiologically complex (e.g. see Fig. 1). It comprises the main intra-nasal passage, the mucous membrane, the ciliary hair-like cells, the mucosal drainage fluid circulating along the internal walls, and the adjoining sinus cavities of various shapes and sizes1,2. Occlusion of the sinus chambers with mucus is associated with many nasal ailments, such as chronic rhinosinusitis3. While surgical treatments essentially focus on enlarging the opening to the sinus chambers, such procedures can be cost-prohibitive, have associated risks, and are mostly reserved for medically refractory diseases. As a first line of treatment, physicians often recommend nasal sprays4,5, with the rationale that these topical drugs will reduce inflammation at the diseased sites and assist in resolving the occlusion and re-establishing natural drainage. However, while such sprays do rank amongst the most commonly used therapeutics, the efficacy of the drugs can be highly patient-specific and there is no well-defined protocol to ensure that precise dosage would reach the intended intra-nasal target sites.

Figure 1
figure 1

(a) Anatomic features inside a human nose, as viewed on an invasive cut-away. (b) Position of the cut-away section, marked by dark line. (c) Representative coronal section, with the main nasal passage shown in lighter color.

Transport of topical drugs inside the nose encounters a number of challenges, namely the airway tortuosity, the sweeping effect of mucociliary drainage, and a lack of consistent usage protocol for the medical devices employed in drug application (primarily owing to the inter-subject heterogeneity in internal anatomic geometry). While optimizing the trajectories for topical nasal drug transport by experimental trials in real human subjects is improbable; with advancement of computational tools, there has been a significant push to obtain numerically simulated predictions of respiratory flow physics and transport therein; see e.g.6,7,8. Of interest are nasal spray simulation studies on in silico models, re-constructed from medical imaging, to measure drug delivery along the nasal passages9, in the sinuses10,11, and on the effects of surgical alterations of the anatomy on nasal airflow12,13,14,15 as well as on topical transport of drugs16,17,18,19. The latter addresses the role of airway channel’s shape in the context of airflow-droplet interactions. Notably, while using medical devices like sprayers, which are inserted at the nostril, the anterior airway geometry gets altered. To simplify the situation though, computational results10 suggest that such initial perturbations do not greatly change or adversely affect the eventual drug deposits at the diseased sites.

Despite the abundance of computational research on nasal drug delivery, there is a distinct lack of articulate instructions for guidance on what could be the “best” way to use the commercially available sprayers. First, numerical studies often do not use a realistic distribution of droplet sizes while simulating topical sprays. Focusing on specific droplet diameters is resourceful while studying the detailed nuances of transport characteristics in that size range; however this somewhat limits the applicability of the subsequent findings while predicting the performance of real sprays, which have a wide variability of droplet sizes in each spray shot. Secondly, the inter-subject anatomic variations also render it difficult to identify a generic spray orientation that can work for all and ensures maximal delivery of drugs at the diseased locations inside the nose.

In this study, we have numerically tracked the transport of therapeutic particulates from over-the-counter nasal sprays via inhaled airflow. The computational fluid dynamics (CFD) models of droplet transport and the in silico prediction of their deposition sites along the nasal airway walls have been compared with in vitro spray experiments in 3D-printed solid replicas of the same anatomic reconstructions. We have proposed a new strategy of nasal spray usage and the recommendation is supported by a significant improvement in target site particulate deposition (TSPD), when compared to the prevalent spray use techniques. The study also expounds20,21,22 on the potential of CFD as a tool in nasal ailment treatment and subject-specific prognosis, and can contribute to the emergence of non-invasive personalized therapeutics and treatment strategies.

Preliminary results pertaining to this work have featured at the American Physical Society (APS) – Division of Fluid Dynamics Annual Meetings23,24 and at the International Society for Aerosols in Medicine (ISAM) Congress25,26,27.

Methods

Anatomic reconstructions

All methods were performed in accordance with the relevant guidelines and regulations, including use of de-identified computed tomography (CT) data from three pre-surgery chronic rhinosinusitis (CRS) patients - collected under approval from the Institutional Review Board (IRB) at the University of North Carolina at Chapel Hill. We also obtained informed consent for participation in this study (which includes obtaining and use of CT data) from the test subjects. Subject 1 was a 41 year-old Caucasian male (body weight 88.0 kg, body mass index 25.3); subject 2 was a 70 year-old Caucasian male (body weight 67.5 kg, body mass index 24.8); and subject 3 was a 24 year-old Caucasian female (body weight 93.1 kg, body mass index 32.6). Medical-grade CT scans of the subjects’ nasal airways were used to re-construct digital models through thresholding of the image radiodensity, with a delineation range of −1024 to −300 Hounsfield units for airspace10,28, complemented by careful manual editing of the selected pixels for anatomic accuracy. As part of that process, the scanned DICOM (Digital Imaging and Communications in Medicine) files for each subject were imported to the image processing software Mimics v18.0 (Materialise, Plymouth, Michigan). For this study, we subsequently considered each side of the nose in the in silico models as a distinct nasal passage model, while studying the droplet transport properties when the spray nozzle was placed on that side: (a) subject 1’s right side constituted nasal passage model 1 (NPM1) and his left side was nasal passage model 2 (NPM2); (b) subject 2’s left side was nasal passage model 3 (NPM3); and (c) subject 3’s right side was nasal passage model 4 (NPM4) and her left side was nasal passage model 5 (NPM5). Note that subject 2’s right-side anatomy did not exhibit a direct access to the diseased intra-nasal targets from outside of the nostril and was not selected for this study. This had to do with the scope of our study design; for details see the section on target site identification. Also refer to the discussion section for follow-up comments.

To prepare the in silico anatomic models for numerical simulation of the inhaled airflow and the sprayed droplet transport therein, the airway domain was meshed and spatially segregated into minute volume elements. The meshing was implemented by importing the Mimics-output in stereolithography (STL) file format to ICEM-CFD v18 (ANSYS, Inc., Canonsburg, Pennsylvania). Following established protocol10,29, each computational grid comprised approximately 4 million unstructured, graded tetrahedral elements; along with three prism layers of approximately 0.1-mm thickness extruded at the airway-tissue boundary with a height ratio of 1.

Full name

Acronym

NPM

Nasal Passage Model

TSPD

Target Site Particulate Deposition/Delivery

CT

Computed Tomography

CRS

Chronic Rhinosinusitis

OMC

Ostiomeatal Complex

CFD

Computational Fluid Dynamics

DICOM

Digital Imaging and Communications in Medicine

STL

Stereolithography

ROI

Region of Interest

NPD

Nozzle Positioning Device

CU

Current Use (spray usage protocol)

LoS

Line of Sight (spray usage protocol)

COVID-19

Novel Coronavirus Disease 2019

List of acronyms.

Inspiratory airflow and sprayed droplet transport simulations

Laminar steady-state models work as a reasonable approximation while modeling comfortable resting to moderate breathing8,30,31,32. Furthermore, with our simulations focusing on a single cycle of inspiration, steady state flow conditions were adopted as a feasible estimate. Based on the principle of mass conservation (continuity), and assuming that the airflow density stays invariant (incompressibility), we have

$$\nabla \cdot {\bf{u}}=0,$$
(1)

with u representing the velocity field for the inspired air. Conservation of momentum under steady state flow conditions leads to the modified Navier-Stokes equations:

$$\rho ({\bf{u}}\cdot \nabla ){\bf{u}}=-\nabla p+\mu {\nabla }^{2}{\bf{u}}+\rho {\bf{b}}\mathrm{}.$$
(2)

Here ρ = 1.204 kg/m3 represents the density of air, μ = 1.825 × 10−5 kg/m.s is air’s dynamic viscosity, p is the pressure in the airway, and b stands for accelerations induced by different body forces. To simulate the airflow, Eqs. (1) and (2) were numerically solved through a finite volume approach, in the inspiratory direction. The computational scheme on ANSYS Fluent v14.5 employed a segregated solver, with SIMPLEC pressure-velocity coupling and second-order upwind spatial discretization. Solution convergence was obtained by minimizing the flow residuals (viz. mass continuity \(\, \sim {\mathscr{O}}{\mathrm{(10}}^{-2})\), velocity components \(\, \sim {\mathscr{O}}{\mathrm{(10}}^{-4})\)), and through stabilizing the mass flow rate and the static outlet pressure at the nasopharynx of the digital models. A typical simulation convergence run-time with 5000 iterations clocked approximately 10 hours, for 4-processor based parallel computations executed at 4.0 GHz speed.

The numerical solutions implemented the following set of boundary conditions: (1) zero velocity at the airway-tissue interface i.e. the tissue surface lining the sinonasal airspace (commonly called no slip at the walls), along with “trap” boundary conditions for droplets whereby a droplet comes to rest after depositing on the wall; (2) zero pressure at nostril planes, which were the pressure-inlet zones in the simulations, with “escape” boundary condition for droplets that allowed outgoing trajectories to leave the airspace through the nostril openings; and (3) a negative pressure at the nasopharyngeal outlet plane, which was a pressure-outlet zone, also with an “escape” boundary condition for droplets. The negative nasopharyngeal pressure was adjusted to generate inhalation airflow rates with less than 1% variation from subject-specific measurements of resting breathing. The physical recordings were collected with LifeShirt vests33 that tracked chest compression/expansion during breathing, and accordingly quantified the inhalation rates (for additional details, see Table 1).

Table 1 This table incorporates the parameters for measured and simulated inhalation airflow, in the study subjects.

After simulating the airflow, sprayed droplet dynamics were tracked through discrete phase particle transport simulations in the ambient airflow, and the corresponding Lagrangian tracking estimated the localized deposition along the airway walls through numerical integration of the following transport equations34:

$$\frac{d{u}_{d}}{dt}=\frac{18\mu }{{d}^{2}{\rho }_{d}}\frac{{C}_{D}\mathrm{Re}}{24}({\bf{u}}-{{\bf{u}}}_{{\bf{p}}})+{\bf{g}}\left(1-\frac{\rho }{{\rho }_{d}}\right)+{{\bf{F}}}_{{\bf{B}}}.$$
(3)

The parameters here are ud, representing the droplet velocity; along with u as the airflow field velocity, ρ and ρd respectively as the air and droplet densities, g as the gravitational acceleration, FB as any other additional body forces per unit droplet mass (as for example, Saffman lift force that is exerted by a typical flow-shear field on small particulates transverse to the airflow direction), and \(18\mu \,{C}_{D}\,\mathrm{Re}\,({\bf{u}}-{{\bf{u}}}_{{\bf{d}}})/24({d}^{2}{\rho }_{d})\) quantifies the drag force contribution per unit droplet mass. Here, CD is the drag coefficient, d is the droplet diameter, and Re represents the relative Reynolds number.

Mean time step for droplet tracking was in the order of 10−5 sec., with the minimum and maximum limits for the adaptive step-size being \( \sim {\mathscr{O}}({10}^{-10})\) sec. and \( \sim {\mathscr{O}}({10}^{-3})\) sec., respectively. Also note that the solution scheme posits the particulate droplets to be large enough to ignore Brownian motion effects on their dynamics. Post-processing of the simulated data laid out the spatial deposition trends, which were then tallied against in vitro observations.

3D printing and physical experiments

To assess the reliability of numerically predicted topical deposition vis-à-vis physical experiments, 3D-printed anatomic replicas were generated for subject 1’s airway and hence included both NPM1 and NPM2. The posterior parts of the solid models were made from the stereolithography material Watershed (DSM Somos, Elgin, Illinois). Post-digitization, the printing job of the posterior component was sub-contracted to ProtoLabs (Morrisville, North Carolina). Printing of the anterior soft plastic part on a Connex3 3D printer was done by Ola Harrysson’s group at North Carolina State University (at the Edward P Fitts Department of Industrial and Systems Engineering), using polymer inkjetting process on Tangogray FLX950 material. See Fig. 2(a–c) for representative pictures of a digitized model and the corresponding 3D replica.

Figure 2
figure 2

(a) In silico model: CT-based digital reconstruction of subject 1’s airway. Panels (b) and (c) respectively show the sagittal and coronal views of the solid 3D-printed replica of the digital model. Note that the solid models comprise a soft outer nose (to mimic the pliability of a real nose) and a posterior hard plastic part. The anterior and posterior 3D-printed components in each model were designed to fit snugly together. Panels (d) and (e) depict the experimental setup for in vitro measurement of sprayed deposits in anatomic solid models.

Recording deposits through gamma scintigraphy

Intra-nasal topical delivery was tracked through in vitro examination of mildly radioactive spray deposits in the 3D-printed anatomic replicas. To ensure that the spray axis orientation and nozzle location aligned with the corresponding simulated spray parameters, we used specially designed nozzle positioning devices (NPD) inserted at the nostril. The spray bottle was fitted into the NPD, while administering the spray via hand-actuation. For each sample test, a bottle of commercial nasal spray Nasacort was labeled with a small amount of radioactive Technetium (Tc99m) in saline. At the time of dispensing the spray shots, a vacuum line controlled by a flow-valve was used to set up inhalation airflow through the model, and the flow rate was commensurate with the subject-specific breathing data (Table 1). Corresponding setup is in Fig. 2(d,e). Four independent replicate runs of each spray experiment were conducted, followed by compilation of the means and standard deviations of the drug deposits along the inner walls of the solid models. The topical deposition was proportional to the radioactive signals emitted from the spray solution traces that deposited inside a solid model and was quantifiable through image-processing of the scintigraphy visuals, collected using a BodyScan (MieAmerica, Forest Hills, IL) 400-mm width by 610-mm height 2D gamma camera. The pixel domain was 256 × 256, with an image acquisition time of 3 minutes; and one pixel equated to a Cartesian distance of 2.38 mm in the digital and 3D models.

Model segmentation for comparison with numerical data

To facilitate the comparison between the numerical predictions on droplet deposition and the physical observation of gamma scintigraphy signals in the corresponding solid replica, we segregated NPM1 and NPM2 into virtual segments oriented along three different directions. Figure 3 lays out the Cartesian coordinate directions for the 3D space. X was perpendicular to the sagittal plane traversing from left to right sides of the nasal models (with the model head facing forward), Y was perpendicular to the axial plane traversing from inferior to superior aspects of the models, and Z was perpendicular to the coronal plane traversing from anterior to posterior aspects of the models. The virtual segments were oriented along the XY (coronal), YZ (sagittal), and ZX (axial) planes. Parallel to the XY coronal plane, the models contained 12 segments (named, C12–C1 sagittal columns); there were 9 compartments (C1–C9 frontal columns) parallel to the YZ sagittal plane, and there were 12 compartments (R1–R12 sagittal rows) parallel to the ZX axial plane (see Fig. 3).

Figure 3
figure 3

Panels (a), (b), and (c) depict the gridline schematic on NPM1 and NPM2, that is used to extract the deposition fractions from the gamma scintigraphy-based quantification of the sprayed deposits in the solid replicas. The models are respectively segregated into 3 sets of compartments: sagittal columns, frontal columns, and sagittal rows. Panel (d) shows the perturbation of the base gridline by 1 pixel. Representative Technetium signals are in panel (e). Note: In regard to the axis system, the circle with solid dot implies out-of-plane direction from this page, the circle with cross signifies into-the-plane of this page.

For each compartment, the particulate deposition fraction predicted from the simulation was compared with the deposition fraction measured based on gamma signals of the deposited particulates in the corresponding compartment of the 3D-printed model. To achieve this, signals emitted from the solution traces, that settled along the airway walls, were subjected to image processing analysis. Therein, by superimposing the compartmental grid on the radio-images, the signals were extracted from each compartment. In order to align the grid on the image in a manner consistent with the virtual model, three inset discs were designed as reference points on the outer surface of the virtual and 3D-printed models. Americium sources from commercial in-home smoke detectors were inserted into the insets as reference points on the 3D-model and a radio-image was recorded. For the analysis, the scintigraphy images were processed using ImageJ35 by constructing a region of interest (ROI) referenced to the fixed Americium sources. Care was taken to align the emitted visual signals with similar reference regions within the superimposed grid. This was done via manual visualization to achieve a best fit of signal intensity within reference regions. The grid compartment planes positioned using this visual best-fit technique were designated as “reference planes”. Given the nature of the radioactive signals and the resolution of the radio-image, some signal intensity resided outside of reference regions even while using best-fit practices. A reasonable fit could be obtained by shifting the image by one pixel in either direction (positive shift/negative shift). In order to account for this variation, alternative plane positions (see Fig. 3(d)) were created by shifting the reference planes one pixel along the positive and negative axes for each set of Cartesian planes. These three sets of compartment planes were positioned in the in silico modeling software using the measured distances from the reference regions. The corresponding Cartesian coordinates of these planes were used to assign droplet deposition locations from the computational simulations to grid compartments, for comparison with the in vitro model. In these comparisons, we left out the deposits in the anterior nose (from the CFD data as well as the physical recordings) in order to negate the bright radiation signal coming from that zone in the experimental deposits; and focused only on measurements from the posterior parts of the respective models. Note that the anterior nose in an in silico model is in fact the removable soft pliable anterior part in the corresponding 3D print (e.g. see Fig. 2).

Identification of target site and spray parameters

Effect of airflow on droplet trajectories

Inertial motion of a droplet is linearly proportional to its mass, and hence is exponentially proportional to the droplet diameter. Consequently, for bigger droplets, the inertial motion persists longer before being taken over by the ambient airflow. Figure 4(a) tracks the trajectory of a representative 5 μ droplet. In there, the tiny red circle marks the location where the inertial motion of the droplet got overwhelmed by the ambient flow, beyond which the droplet trajectory was same as the airflow streamline on which it was embedded at the red circle’s location. Note the contrasting 25 μ droplet trajectory in Fig. 4(b), where the inertial motion persisted longer. The phenomenon has a significant impact on drug deposition trends. The bigger droplets (≥100 μ) show a greater propensity to hit the anterior walls directly owing to their high initial momentum, while smaller droplet sizes penetrate further into the airspace; see e.g. Figure 4(c,d). To ensure that the bigger droplets also reach the target sites, we argue that it would be conducive to harness their inertial motion and direct those droplets actively toward the target when they exit the spray nozzle. This can be feasibly achieved by orienting the spray axis to pass directly through an intended anatomic target zone.

Figure 4
figure 4

Comparison of representative trajectories for a 5 μ droplet and a 25 μ droplet in a sample sinonasal airspace. In panel (a), the smaller droplet has weaker inertial momentum and the ambient airflow streamline takes over its motion much earlier than that in case of a heavier droplet like the one in panel (b), where the inertial momentum of the 25 μ droplet persists longer. The small red circle in (a) depicts the point where the inertial momentum gets overwhelmed by the fluid streamline. Evidently, owing to smaller inertia, the droplets with smaller diameters get predominated by the airflow streamlines earlier than the bigger droplets. This results in a better penetration and spread of sprayed droplets in the nasal airspace, as shown in panel (c), for a different nasal model. On the contrary, spray shots with exclusive share of bigger droplets (e.g. ≥100 μ here) tend to follow their initial inertial trajectories, without much effect of the airflow streamlines on their paths, and deposit along the anterior walls of the nasal airspace, as depicted in panel (d). The red boundaries in panels (c) and (d) highlight the difference in particulate penetration into the model, in the two cases. Note: These images were created using FieldView, as provided by Intelligent Light through its University Partners Program.

Current use instructions

Inconsistency and ambiguity in instructions36,37 indicate a lack of definitive knowledge on the best ways to use a nasal spray device. Different commercial sprayers often offer somewhat contrasting recommendations. However, there is a common agreement (see Fig. 5(a)) that the patient should incline her/his head slightly forward, while keeping the spray bottle upright36,38. Furthermore, there is a clinical recommendation to avoid pointing the spray directly at the septum (the separating cartilaginous wall between the two sides of the nose). These suggestions were adopted in our standardization16 of “Current Use” (CU) protocol for topical sprays. The digital models were inclined forward by an angle of 22.5°, and the vertically upright36 spray axis was closer to the lateral nasal wall, at one-third of the distance between the lateral side and septal wall. Also, the spray bottle was so placed that it penetrated into the airspace by a distance of 5-mm, inspired by the package recommendations of commercial sprayers38 for a “shallow” insertion into the nose. Refer to Fig. 5(b,c) for the schematics of the CU protocol used in this study.

Figure 5
figure 5

(a) Illustration depicting nasal spray usage, as per instructions38 that come with over-the-counter sprays. Panel (b) and inset (c) depict the protocol implemented in the numerical simulations for the “Current Use” (CU) spray orientation. Note that δ is the linear distance between lateral wall and septum (the cartilaginous “mid-wall” in the nose, separating right and left airways) at 5-mm insertion depth into the nose. The model “head” is tilted forward by 22.5°. The vertically upright dashed line represents the spray nozzle axis.

Target site identification and proposing an alternate spray use criteria

All sinuses, except sphenoid, drain into the ostiomeatal complex (OMC), it being the main mucociliary drainage pathway and airflow exchange corridor between the nasal airway and the adjoining sinus cavities. To ensure that as many drug particulates reach the sinus chambers and their vicinity as would be possible, we hypothesize that the spray axis should be directed straight toward the OMC39. This is supported by our observation of the effect of airflow physics on droplet trajectories. If the spray axis hits the OMC directly, the likelihood that the larger droplets will deposit there is higher. We refer to this usage protocol as “Line of Sight” (LoS). Like the CU protocol, the LoS protocol also had the sprayer inserted at a depth of 5-mm into the nasal airspace. Representative LoS orientation is shown in Fig. 6.

Figure 6
figure 6

Panels (a) and (b) show the locations of the main target sites in a representative sinonasal reconstruction, i.e. the OMC (acting as the mucociliary drainage pathway for the sinuses) and the sinus cavities. Panels (c)-(e) demonstrate the “Line of Sight” (LoS; represented by the black lines) in NPM1. The anatomic zone, colored red, marks the OMC. Note that panel (d) is the 3D-printed soft nose from NPM1, exhibiting the same approximate orientation as that of the digital model in panel (c), giving a direct straight-line access to the target sites, and hence an LoS. The blue component in the image on panel (d) indicates the approximate location of the OMC.

TSPD percentage at the OMC and the sinuses was evaluated as = 100 × (Mtarget/Mspray); with Mtarget being the spray mass of the particulate droplets deposited at the OMC and inside the sinus cavities, and Mspray being the mass of one spray shot.

Generation of varying peripheral directions around the true CU and LoS directions

To establish the robustness of the TSPD predictions for the CU and LoS protocols, we also tracked droplet transport and deposition when the spray directions were slightly perturbed. Such perturbed peripheral directions for CU initiated 1 mm away on the nostril plane and were parallel to the CU’s vertically upright true direction. For LoS, the perturbed peripheral directions were obtained by connecting the base of the true LoS direction on the nostril plane with points that radially lie 1 mm away from a point on the LoS; this specific point being 10 mm away along the LoS from the base of the LoS direction on the nostril plane (e.g. see bottom panel of Fig. 7 for an illustrative example).

Figure 7
figure 7

Comparison of the simulated spray deposits from the CU and LoS protocols. The yellow bars represent the TSPD for the CU spray orientations, and the blue bars quantify the TSPD recorded for the LoS spray orientations. The gray bars are the predicted deposits when the true CU and LoS directions were perturbed by 1 mm. Panels (a)–(e) are the results for five different airway models: Nasal Passage Model 1 (NPM1), Nasal Passage Model 2 (NPM2), Nasal Passage Model 3 (NPM3), Nasal Passage Model 4 (NPM4), and Nasal Passage Model 5 (NPM5). Panel (f) compares the TSPD for peripheral directions in a 0.5-mm perturbation (on the left) with respect to a 1-mm perturbation (on the right) from the true LoS orientation, both in NPM1. As expected from the overall findings, the TSPD increased for the perturbed spray directions that were closer to the true LoS. Panel (g) depicts the spatial perturbation parameters for the LoS spray axis orientation in NPM1.

Parameters for the simulated spray shot

Over-the-counter Nasacort (Triamcinolone Acetonide), a commonly prescribed and commercially available nasal spray, was selected for this study. Four units of Nasacort were tested at Next Breath, LLC (Baltimore, MD, USA) to characterize the in vitro spray performance. Corresponding plume geometry was analysed through a SprayVIEW NOSP, which is a non-impaction laser sheet-based instrument. Averaged spray half-cone angle was estimated at 27.93°, and the droplet sizes in a spray shot followed a log-normal distribution. With the droplet diameter as x, the droplet size distribution can be framed as a probability density function of the form40:

$$f(x)=\frac{1}{\sqrt{2\pi }x\,\mathrm{ln}\,{\sigma }_{g}}\exp \left[-\frac{{(\mathrm{ln}x-\mathrm{ln}{x}_{50})}^{2}}{2{(\mathrm{ln}{\sigma }_{g})}^{2}}\right].$$
(4)

Here, x50 = 43.81μ is the mass median diameter (alternatively, the geometric mean diameter41) and σg = 1.994 is the geometric standard deviation. The latter quantifies the span of the droplet size data. Measurements were also made with and without the saline additive in the sprayer, and the tests returned similar droplet size distribution. Note that a saline additive was used during the physical recording of the sprayed deposits. Also, as per earlier findings in literature42, the mean spray exit velocity from the nozzle approximates at 18.5 m/s, based on phase doppler anemometry-based measurements.

For the test spray units (at Next Breath), the actuation forces were found to range between 8.5–11.0 kg-force. Considering an actuation area of approximately 4 cm2, the force measurements agree well with earlier values in literature43,44,45 and hence the resultant pressure exerted on the droplets in our physical experiments was assumed to maintain a similar droplet size distribution, as was determined in the test cases by Next Breath. The droplets contained in one spray shot in the numerical simulations followed the same size distribution.

While simulating the droplet trajectories, we assumed typical solid-cone injections and tracked the transport for 1-mg spray shot while comparing the TSPD trends from the CFD predictions with the corresponding experimental drug delivery patterns. On the other hand, 95.0306 mg (which is one shot of Nasacort, as quantified by Next Breath, LLC) of spray mass transport was simulated while comparing the CFD-based TSPD numbers for the LoS and CU protocols in each model.

Results

Comparison between CU and LoS spray usage protocols

LoS was found to be consistently superior in comparison to the CU spray placement protocol, while targeting the OMC and the sinus cavities for drug delivery. Table 2 lists the deposition fraction percentages for each spray release condition in the five airway models (NPM1–NPM5). For a graphical interpretation, we have plotted the same information on Fig. 7. Overall, the deposition fraction for the LoS was on an average 8.0-fold higher than the CU deposition fraction, with the corresponding subject-specific improvement range being 1.8–15.8 folds for the five test models. The improvement does decay when the perturbed peripheral spray directions are compared, to assess the robustness of the LoS protocol’s advantage over CU. Considering the varying peripheral directions around the true LoS and CU, the LoS set registered an average 3.0-fold increase in TSPD, with the corresponding subject-specific improvement range being 1.6–4.3 folds.

Table 2 Numerical prediction of targeted drug delivery from CU and LoS protocols.

Statistical tests – on improvements achieved by the revised spray use strategy

LoS was compared to CU through a paired study design on the data from five test models. Table 3 lays out the computed numbers. For each model, the outcome comprised the percentage of deposition in OMC and the sinuses for both CU and LoS spray usage. Null hypothesis considered for this statistical test assumed that the TSPD would be same for CU and LoS in an airway model. The deposition percentage corresponding to CU and LoS protocols in the same nostril were treated as paired observations for a paired t-test to check the null hypothesis. Owing to a relatively small study cohort, paired Wilcoxon signed rank test was also used for robustness check. In order to study how spatial variation might affect the difference between CU and LoS, three different ways of calculating the percentage of deposition were implemented. The first strategy considered the average deposition from the true LoS and CU directions. The second strategy compared the TSPD averaged from the true CU and LoS directions, along with the deposition data for spray release parameters obtained by perturbing the respective true directions. The third strategy used TSPD averaged exclusively from the deposition data corresponding to the perturbed spray release parameters. This allowed us to assess the robustness of any probable improvement from using LoS, while still accounting for slight spatial variations of the spray direction.

Table 3 Statistical tests for the comparison between CU and LoS protocols.

The first comparison method demonstrates an average deposition increase of 5.4 percentage points for LoS (6.39-% for LoS vis-à-vis 0.98% for CU). This difference is significant at the 0.05 level with a p-value from the paired t-test of 0.03. The paired Wilcoxon signed-rank test has a p-value of 0.06, which was the lowest possible p-value for the Wilcoxon signed-rank test given only five pairs of data. In the second comparison scheme, LoS has an increased deposition of 1.62 percentage points relative to CU (2.49% vis-à-vis 0.87%). The p-value for this difference is 0.02 using the paired t-test and 0.06 using the Wilcoxon signed rank test. Finally, for the third comparison method, LoS registered an increased deposition of 1.05 percentage points relative to CU (1.90% vis-à-vis 0.86%). The p-value for this difference is 0.02 using the paired t-test and 0.06 using the Wilcoxon signed rank test. This provides a strong evidence that LoS leads to higher percentage of deposition in the OMC and sinuses. The estimated difference is largest when using just the true directions, but the difference is still statistically significant even when using the spray release points obtained by perturbing the true directions. The p-value from the paired t-test is actually lower when the TSPD from just the perturbed points are considered, owing to the reduced variance for the estimated difference. For all three ways of estimating the percentage of deposition, the paired Wilcoxon signed-rank test returns a p-value of 0.06. With only five pairs of data, this suggests that the use of LoS does result in statistically significant higher deposition for all five nostril models.

Comparison of the simulated TSPD predictions with physical experiments

Figure 8 compares the numerical TSPD predictions with corresponding gamma scintigraphy-based experimental recordings in NPM1 and NPM2. While the compartmental deposits visibly presented a congruous trend in the sagittal columns, sagittal rows, and frontal columns; we conducted additional statistical tests to verify the homogeneity between the two sets of data so as to establish the reliability of the computational findings.

Figure 8
figure 8

(a) Comparison of the numerically simulated compartmental findings in Nasal Passage Model 1, with respect to the gamma scintigraphy recordings from the corresponding 3D-printed replica. (b) Comparison of the numerically simulated compartmental findings in Nasal Passage Model 2, with respect to the gamma scintigraphy recordings from the corresponding 3D-printed replica. The blue “reference” lines trace the CFD predictions for TSPD in each compartment, with the light gray and dark gray lines respectively marking the variability in prediction, for +/−1 pixel shift while superimposing the gridlines on the numerical data-space. The yellow lines trace the TSPD recorded from the physical experiments.

Table 4 gives the Pearson and Kendall’s correlation between the numerical and experimental models for the average deposition fractions in NPM1 and NPM2 for the LoS protocol. The confidence intervals are based on 1000 bootstrap samples, instead of asymptotic approximations, because of the relatively small sample size. Based on the output, we can see that the Pearson correlation is consistently very high while the Kendall’s correlation is somewhat lower. However, while the Kendall’s correlation is frequently thought to be more robust to outliers, particularly for small sample sizes like this data-set; in this particular instance the Pearson correlation is likely more illustrative. This is because the Pearson correlation is able to show that, for the most part, the magnitudes of the estimates are similar and comparable between the numerical and experimental models. In general, there is a strong linear relationship between the percent of deposition prediction from the numerical model and the corresponding physical measurements in the experimental model. The lower Kendall’s correlation (overall mean measure 0.78) is largely due to regions where both the numerical and experimental models had very low average deposition but the exact rank of these regions changed considerably between the two data-sets. Note that this does not necessarily indicate a poor performing numerical model. However, the relatively high Pearson correlation (overall mean measure 0.91) does indicate that the numerical models perform well while predicting the sprayed droplet transport.

Table 4 Comparison between the compartmental data from numerical simulations and physical experiments.

Discussion

CFD-guided nasal spray usage defined by the LoS protocol was found to significantly enhance topical drug delivery at targeted sinonasal sites, when compared to currently used spray administration techniques. With increased sample size, this work can be the catalysis toward prompting personalized instructions and specifications for improved use of topical sprays. The findings, thus, have the potential to substantially upgrade the treatment paradigm for sinonasal ailments through the ability to ascertain LoS in individual subjects via endoscopic examinations conducted in the clinic, and to help guide treatment decision-making and patient instructions for spray usage.

Concept of LoS scoring and on the adaptability of our findings in clinical practice

As means to quantifying the suitability of a person’s airway for the LoS spray protocol, we exploratorily propose a scoring system that is based on how much of the targeted drug delivery sites (OMC, sinuses) are visible when inspected clinically from outside of the nostril. The scoring system will also serve to quantify nasal anatomic variability among individuals. Accordingly, as part of the current study, the LoS scores (see Table 5) were first determined observationally, based on the external visibility of the OMC site in the in silico sinonasal reconstructions. We fixed a range of scores [1, 4], with 4 being used when the LoS direction was easiest to ascertain. Subjective as that scoring procedure may be, it is similar to what attending physicians will gauge during a clinic visit to determine if a particular patient has a “line of sight” in her/his nasal anatomy. So, to establish the relevance of the findings from this manuscript toward revisions of the therapeutic protocol for sinonasal care, it is important to assess the comparability of the observational LoS scores with more objective score determination techniques. This was achieved by calculating the surface area of the nostril plane and the projected area of the OMC on the plane of the nostril. We computed the ratio of the projected area to the nostril area, as a percentage. Scores of 4 were assigned if the ratio exceeded 6%, 3 if the ratio exceeded 4%, 2 if the ratio was more than 1.5%, and 1 if the ratio was greater than 0%. The two scoring techniques yielded very similar results (as in Table 5), with the highest and lowest scores respectively going to the same anatomic models. Pearson’s rank correlation for the two sets of scores was 0.85. While a broader study, involving clinical trials, will be necessary to revise therapeutic protocol for nasal drug delivery, the present results illustrate the easy adaptability of our findings into clinical practice settings.

Table 5 Comparison of the LoS scores, obtained observationally and through determining the surface area projection of the targeted OMC on the nostril plane.

On the comparability of the experimental data with the numerical findings

The computational simulations assumed a laminar framework to mimic steady breathing. However, one may argue that even with resting breathing rates, the airflow often contains transitional features like vortices, emerging from the roll-up of shearing fluid layers during flow-structure interactions46,47,48,49,50,51 at the anatomic bends. Some of these nuances are, in fact, difficult to model without proper turbulence simulations52,53. However, true as that may be, the effect of these flow artifacts on eventual drug delivery in the sinuses has been found to be somewhat nominal while comparing laminar and turbulence simulation results10.

On the other hand, the in vitro techniques also often pose challenges. For instance, there can be post-deposition run-off as the deposited solution traces undergo translocation along the inner walls of the solid replica. Such drip-off dynamics can lead to a flawed estimate of regional deposition. The effect of post-deposition dripping can be conjectured to be most prominent for the signals extracted from the sagittal rows, as the deposited droplets start moving downward along the internal solid walls of the 3D-printed models, owing to gravitational effects. This is confirmed by the physical and numerically-predicted signals from the sagittal rows demonstrating relatively lower correlation coefficients (when contrasted with the correlations for the signals from the sagittal and front columns) in the two experimental comparisons (e.g. see Table 4).

In the gamma scintigraphy-based method of recording deposits, the radiation signal undergoes some level of scattering and hence in the process of signal extraction from each of the compartments, there is the possibility that signals from one compartment may contaminate the signals at neighboring compartments. To minimize this effect while carrying out the comparisons, the nose (the soft plastic anterior part in the 3D-printed models), which had a bright radiation signal owing to the relatively large amount of anterior deposits, was excluded from both the experimental and numerical data.

Finally, while the inhalation airflow rates were same in vitro and in silico, the airflow partitioning on the two sides of the nasal airways was likely affected by the placement of the NPD, while administering the spray through hand-actuation.

Caveats and future implications

Readers should note that this was a computational study with validation from spray transport observations in inanimate solid replicas. Also, not every patient will have a clear access to the OMC, and hence may be without an LoS. For instance, in the current study, of the six airway sides in the three study subjects, subject 2’s right-side airway did not exhibit an LoS.

Bulk rheology of the spray also affects the droplet size distribution. The spray property measurement tests having been performed in real over-the-counter sprays, we did not separately examine the effect of different droplet viscosities on the spray deposition trends. A different viscosity of the nasal spray can indeed alter the drug deposits, as observed in multiple studies37,54,55. It should however be pointed out that the spray positioning strategies proposed in this study could be conjectured to be generic and should maximize drug delivery to the OMC and the sinuses for other sprays as well.

It is also critical to note that the flow simulations for evaluating the spray usage strategies were not multiphase; implying that the sprayed droplets were not affected by constituents such as inhaled air moisture and the mucous lining, nor was there any consideration of droplet evaporation. There are, however, earlier findings in literature that have looked at some of these nuances; e.g., on the interaction of deposited particulates with mucus56 and on the phase change of inhaled droplets during their passage through the respiratory tract57. The current study simply tracked the motion of inert droplets against the ambient inspiratory airflow and recorded their regional deposition. Based on the surfactants in the spray solution, the droplets might also be rendered hydrophilic; but such effects are beyond the scope of this project and the numerical schemes that have been implemented. It is, however, expository to reckon that such hydrophilicity may at times lead to agglomeration of droplet molecules, which can impact the topical drug deposition estimates.

This study, its restricted sample size and limitations notwithstanding, is still, to the best of our knowledge, the first-of-its-kind to propose an alternative easy-to-implement strategy that can significantly improve the intra-nasal delivery of topical drugs at the diseased sites. The recommendation for using the “line of sight” is user-friendly, personalized (the physician can instruct the patient on the spray usage technique based on a fast LoS check in the clinic), and has the potential to be smoothly incorporated into the nasal standard-of-care. For probable revisions to the clinical regimen, we will need a broader study with more subjects, along with a component for clinical trials to track patient response. Comparison of the numerical data with in vivo spray performance will also eliminate errors that contaminate the in vitro TSPD numbers (e.g. from drip-off of the deposited solution along the inner wall contours of the 3D-printed models). Nevertheless on a larger intriguing perspective, the current study conclusively postulates how relatively simple engineering analysis and mechanistic tools can usher in transformative changes in the prognosis and treatment protocol for ailments such as nasal congestion and respiratory infection.

Special comments on the significance of the findings in view of the 2019–2020 Coronavirus pandemic

With the rapid spread of the Novel Coronavirus Disease 2019 (COVID-19) worldwide, it is essential that a vaccine or a curative is developed at the earliest. With respiratory mucosa as the initial site in coronavirus infection and transmission; mucosal immunization through targeted intra-nasal vaccine promises to be an effective strategy for prophylaxis, by inducing mucosal and systemic immune responses. As of May 2020, several research groups are working on the possibility of designing intra-nasal vaccines for COVID-1958,59,60, with supporting data from work carried out on earlier strains of coronavirus61. In this context, the intra-nasal anatomic targeting strategies (e.g. see Fig. 9) discussed in the current study can be of significant help to increase the topical delivery.

Figure 9
figure 9

(a) Superimposed images portraying a representative sagittal comparison of the newly-devised LoS direction, contrasted with the spray bottle positioning (based on instructions that come with over-the-counter nasal sprays) by a volunteer subject. The trial was non-invasive, and the volunteering subject was one of the co-authors of this study. (b) Endoscopic view of the volunteer’s internal nasal anatomy when inspected along the LoS direction. These graphics can assist the readers in adopting the new spray usage strategies with relative ease, and can be of help to improve the application of nasal drugs; especially in the context of developing intra-nasal vaccines for COVID-19.