Use of Magnetic Nanoparticles as In Situ Mucus Property Probe

: Magnetic nanoparticles (MNPs) are unique in their abilities to penetrate and interact with a wide range of liquid media. Because of their magnetic properties, MNPs can be directed to any area of interest, and interact with core structures deep inside the medium which is normally inaccessible. In this report, we investigate the behavior of MNPs in a speciﬁc biological ﬂuid, namely in a mucus layer of air–liquid interface cultured primary normal human tracheobronchial epithelial cells. Using Fokker–Planck algorithm simulations and observing the behavior of MNPs from prior experiments, we found MNPs that are initially less than 100 nm in size, to aggregate into sizes of ∼ 50 µ m and to deviate from the expected Fokker–Planck distribution due to the mucus structure. Based on our analysis, human tracheobronchial epithelial (NHTE) cell mucus viscosity ranges from 15 Pa · s to 150 Pa · s. The results not only conﬁrm the possible use of MNPs as a means for medical drug delivery but also underline important consequences of MNP surface modiﬁcations.


Introduction
Nanoparticles (NPs), both natural and synthetic, have been used as devices for material modification and micron sized molecule transport in a wide range of applications. Fields of application include material science [1], vaccines [2], and its use as vehicles that travel inside living cells and organisms on existing biological transport systems [3,4]. Modifications on the surface of NPs allow researchers to utilize chemical interactions and carry out various tasks. Such interactions between NPs and biological fluids or among NPs themselves reveal details of the biological fluid's microenvironment within the vicinity of each NP. As researchers continue to develop new ways of modifying and using NPs, the potential for their application across physical and biological fields is vast.
A closer look at NP behavior in various biological fluids might reveal how they interact physically and help enhance future applications. It is especially important to understand the underlying physical laws that govern NP parameters such as size and travel speed [5].
In general, small particles of less than 100 nm in size moving in viscous fluids experience forces characterized by small Reynolds numbers [6]. In such environments, drag force is linearly dependent on velocity, and the flow is laminar. In viscous fluids, Stokes' law can be applied in situations where the flow is laminar, particles are smooth, spherical, and do not interact with one another [7]. However, when particles do interact, Stokes' law has to be modified, such as in the case of ions moving through aqueous solutions, as discussed in Ref. [8]. In this study, an imagined electric field is placed around an ion of certain charge, Ze. The ion is imagined to be either stationary or moving at a constant speed in a fluid, so that the electrical force balances drag force. This equation then allows canceling of ion speed and expresses the ion radius in terms of Ze and fluid viscosity, η. The study found discrepancy between this calculated radius and experimentally observed radius. In practice, ions moving through fluids may drag other molecules along with them as they move, due to ionic interactions. The result of this is an effective particle radius, called Stokes' radius, that arises from ionic mobility effects [9] but is distinct from the ion's physical radius. This Stokes radius is a direct result of modified Stokes' law, taking into account ionic charges and their effect on surface adhesion. Another example is marine plankton aggregates [10]. Phytoplankton, the principal source of energy for marine life, form aggregates with minerals and zooplankton when sinking in the ocean. Because minerals make phytoplankton surfaces round and smooth, phytoplankton-mineral aggregates follow Stokes' law, while plankton-only aggregates follow the full Navier-Stokes equation, or else have to be simulated with a modified Stokes' law. Modification to Stokes' law takes into account the molecular aggregate's porosity in the case of marine aggregates, or the degree of adhesion at the particle surface in the case of ionic particles. Often, these modifications reflect some aspects of the particle or fluid's chemical or physical properties, thus we might perceive deviations from Stokes' law as an indication of unrevealed interactions between particles.
Mobile NPs also reveal other useful physical parameters of the biological fluid. In highly viscous environments, a strong drag coefficient forces a moving particle to travel at its terminal velocity almost immediately. The time that it takes to reach this terminal velocity is easily calculated [11]. This characteristic time is very small (∼ 10 −6 s) in viscous environments, and thus terminal velocity can be readily observed in experiments. With this observed terminal velocity, the viscosity of the fluid can be estimated. Other parameters of the fluid, such as pH, density, and surface tension, can also possibly be probed by NPs, using other appropriate physical laws, such as the Navier-Stokes equation [12].
Magnetic NPs (MNPs) with various surface modifications are especially suited for medical purposes such as targeted drug delivery and magnetic resonance imaging [13][14][15]. Magnetic properties of these NPs add an extra level of control for traveling speed, direction, and particle size. Calculation of MNP aggregate size is helpful in experiments that are limited by the physical apparatus. Optical microscopes are limited by the diffraction limit of visible light, which sets the resolution to approximately 200 nm. Aggregates smaller than this size cannot be resolved. Moreover, in physical size determination techniques, such as dynamic light scattering (DLS), NPs need to be well dispersed in a solvent and stable overtime to produce the Brownian motions necessary for characterizing NP sizes from scattered light. With MNPs that have motions directed by magnets and positions described by time-dependent probability density, it becomes necessary to remove MNPs from the biological fluids and into a sample to assess their sizes by DLS. However, doing so not only alters their chemical environment but also removes the physical magnetic field, both of which greatly affect aggregate sizes [16]. In the presence of an applied magnetic field, magnetic dipoles are induced causing MNPs to aggregate. Therefore, removing this induced magnetization is especially detrimental to aggregate size determination. In the present study, assuming MNPs have no ionic charges that interact with their surrounding fluid, we were able to estimate aggregate sizes with Fokker-Planck simulations.
Here, we use the Fokker-Planck equation to simulate MNP movement across corn syrup and the mucus layer of air-liquid interface cultured primary normal human tracheobronchial epithelial cells under a controlled magnetic field [17]. The Fokker-Planck equation includes velocity terms and a diffusion term [18]. The diffusion term involves Stokes' law, and the velocity terms are due to gravity and a magnetic force. From the simulation of the diffusion term, deviation from Stokes' law in a prior research [17] is detected. From the observation of particle distributions, we assume that Stokes' law is valid, but the resulting distribution is distorted by the mucus structure. Gravity and magnetic force are similarly calculated and simulated in the velocity terms of the equation. Comparing these simulations to previous studies [17], we notice that the magnetic force is the major driving force for traveling speed. Therefore, in a high gradient magnetic field, the MNPs will have a chance to penetrate mucus before rapid mucus clearance removes them [17]. We also notice that the MNPs were falling at a much faster speed than that calculated from gravity and the magnetic field. Since there are no other external forces present, we estimate that the radii of the particles were much larger than their initial size. Initially about 20 nm in diameter, MNP sizes ranged from 10-100 microns in diameter, depending on the solvent (biological fluid), type of MNPs, and presence of an external magnetic field. Therefore, we conclude that MNP aggregates similar to those described by [19] had occurred in these experiments. The presented studies could be performed in situ, which is the preferred mode of observation whenever physical and chemical fluid environments greatly affect NP behavior.

Characteristic Time and Terminal Velocity
Under gravity, the time it takes to reach terminal velocities (characteristic time) was estimated by solving the following differential equation [11]: where γv is the drag force, and γ is the drag coefficient from Stokes' law. mg is the gravitational force with m being the mass of the particle and g denoting the gravitational constant near Earth's surface. The solution v(t) is an exponentially converging function that reaches the terminal velocity, v ter , at large times, with an exponential time constant on the order of m/γ. This is called the characteristic time, τ, given as Using the same principles under a magnetic field, Equation (1) becomes where M sat V dB(z) dz is the magnetic force. M sat represents magnetic saturation of MNP, V its volume, and dB(z) dz the gradient of the magnetic field B(z) with respect to z. Here, z measures the distance from the magnetic tip used to pull the MNPs through the fluid [17], see Figure 1 for the experimental setup. The solution v(t) is again exponentially converging to v ter , with the same exponential time constant τ. In the case of gravity, it is easy to see that v ter is reached within a few τ, and the NP sinks at a constant speed afterward [11]. It is not so obvious in the case of a high gradient magnetic field. Here, the magnetic force is a function of position due to the B(z) term.
The gradient, dB(z) dz , increases as the particle sinks closer to the magnetic tip, so v ter would increase as the particle sinks. The particle would only reach v ter if its speed is fast enough so that v ter does not increase significantly over the time τ. Simple calculations show that this is true, and that it is considered acceptable to neglect the small τ. This same assumption was also used in a similar research [20].
Since τ of a MNP is negligible, v ter can be observed directly in certain experiments. This can be estimated in one of two ways. First, using data from a previous study [17], the approximate time for MNPs to fall through the specified thickness of corn syrup or mucus is used as a rough estimate of its v ter . Second, in the present study, MNPs were allowed to sink under both gravity and magnetic field in corn syrup. Then, a plot of MNP position vs. time was used with appropriate curve fit to more accurately determine v ter of MNPs graphically. All Fokker-Planck simulations were based on the first approach, but the second graphical approach can be used to reveal slight differences in the behavior of different MNP types.

Aggregate Radius
Terminal velocity can also be calculated analytically by balancing the drag force with all external forces (setting a = 0 in Equation (3)) [21]: Here, as before, mg is the gravitational force, M sat V dB dz is the magnetic force, and γv ter the drag force. The drag coefficient γ is defined by Stokes' law: where F d is the drag force, v ter is the terminal velocity, and γ is the drag coefficient defined by γ = 6πRη.
Here, η is the fluid viscosity, and the total radius R is the sum of magnetic aggregate radius R m and non-magnetic particles aggregate radius R n : Non-magnetic particles (NMPs) are NPs that are not pulled by an applied static magnetic field (MF) and have negligible mass relative to MNPs. Both gravity and MF have little effect on NMPs, but they can form aggregates of considerable size causing drag force to increase significantly. NMPs include surface modification molecules such as AP3, PLGA, and Dex, etc.
To estimate R, we first calculated its v ter with only one MNP per aggregate using Equation (4), or The discrepancy between experimental and calculated terminal velocity was observed, and the number of MNPs (N m ) and NMPs (N n ) in one aggregate was increased in a 100 to 1 ratio until their v ter s were in agreement. The ratio of 100 magnetic unit to 1 non-magnetic unit was estimated from observations of TEM images in Ref. [17]. From these numbers and the radii of single MNPs (r m ) and NMPs (r n ), the total radius R was calculated by Equation (7) using and Of special importance in force calculations was the increased magnetic radius, R m . It represented the MNP's magnetic aggregate radius, where many tiny MNPs with radii r m are assumed to clump together to form the magnetic part of the aggregate. Thus, the volume V and mass m of the gravitational and magnetic forces in Equation (4) can be obtained by and where ρ is density of the MNPs.

Terminal Velocities: Graphical Approach
The terminal velocity data points were found experimentally using two different methods. In one method, the v ter s were determined with the presence of gravity only. In a second method, the v ter s were determined with the presence of both gravity and a magnetic force pulling down as well.
In the first method, the particle clusters were pipetted into the cuvettes of corn syrup. Sufficient time was taken to ensure that the particle clusters had settled and were no longer impacted by the flow in the corn syrup created by the pressure of the pipette. The particle clusters were pipetted in such a way that they were below the surface tension on the corn syrup. Once it was determined the particle clusters were only moving based on the impact of gravity alone, their position, in mm above from the bottom of the cuvette, was recorded. An individual particle cluster that showed distinctness from other clusters and was spaced far enough away from other clusters, making it easier to track, was selected. A stopwatch was used to record time. The position of the particle cluster was recorded every 5 min until the cluster reached the bottom of the cuvette, or until 50 min had passed. This experiment was done with three different cuvette samples for each type of particle (BaFe 12 O 19 NPs and Fe 3 O 4 NPs), and terminal velocity was averaged for each.
In the second method, the particle clusters were pipetted into the cuvettes of corn syrup similar to how they were in the first method. Once again, sufficient time was taken to ensure that the particle cluster had settled and were no longer impacted by the flow in the corn syrup created by the pressure of the pipette. The particle solution was pipetted into the corn syrup, below the surface just as they were in the first method. At this point, a particle cluster was selected, and its position was recorded. For each trial, the cuvettes were then placed in a fitted square dugout, on top of the apparatus, as seen in Figure 1, behind the z-axis line. A stopwatch was used to record time. The cuvette was removed from the apparatus every 5 s, and the new position of the cluster was recorded. Attention was paid, to ensure the same cluster was being recorded for position each time. Just as was done during the first method, this experiment was done with three different cuvette samples for each type of particle (BaFe 12 O 19 NPs and Fe 3 O 4 NPs), and the measured terminal velocity was averaged for each.
After the experiment, position vs. time data points were plotted with the appropriate curve fit, and v ter was estimated from the plot. For each MNPs, three trials to estimate v ter were done, and the average is taken with uncertainty being half of the difference between the highest and the lowest velocity. One of these trials for gravity and MF is shown in Figures 2 and 3, respectively.

Magnetic Nanoparticle Types and Surface Modifications
In the previous study [17], NPs derived from iron oxide, Fe 3 O 4 , nanoparticles (FeNPs) and barium hexaferrite, BaFe 12 O 19 , and nanoparticles (BaNPs) were allowed to fall through corn syrup and normal human tracheobronchial epithelial (NHTE) cell mucus under gravity and high gradient magnetic fields. FeNps and BaNps are widely used as biological MNPs because of their unique properties. FeNps are superparamagnetic, which makes them react very strongly to MFs, especially when they form clusters [14]. Some FeNps were also selected for their good biocompatibility and low toxicity [13]. BaNps have similar magnetic properties to FeNps, but BaNps are five times more isotropic. This makes them react very strongly to alternating (AC) magnetic fields, making BaNps ideal for cancer therapy techniques such as AC fluid hyperthermia [17,22]. Confocal laser scanning microscopy showed MNP concentration as a function of time as the MNPs traveled through mucus under gravity and magnetic field. The position of the maximum of this concentration represents the average position of the MNPs, see Figure 4. The speed at which this maximum moves is taken as the terminal velocity in mucus. Furthermore, each MNP type was surfaced modified for dispersion and drug delivery purposes. BaNPs, with initial radius of 0.250 µm, were milled to decrease its size to ∼0.010 µm and modified with DL-2-amino-3-phosphonopropionic acid (AP3) to allow bonding to poly(lactic-co-glycolic) acid-dexamethasone (PLGA-Dex). PLGA acts as a carrier for the anti-inflammatory drug, dexamethasone. Further modification with chitosan FITC allows fluorescent detection. FeNPs or polyvinylpyrrolidone coated FeNPs (PVP-FeNPs), with initial radii 0.015 µm and 0.025 µm respectively, were similarly modified with PLGA-Dex-AP3-FITC (without milling), and subjected to the same magnetic tests. Initial sizes of Mnps and other molecules are listed in Table 1. To see the effects of surface modifications on different types of MNPs in this study, we simulated the behavior of each surface-modified BaNP and FeNP separately. Table 1. Radii of single magnetic nanoparticles (MNPs) and surface modification molecules (r m and r n ). Particle sizes and uncertainties were estimated from transmission electron microscopy (TEM) images or size distribution plots given in Ref. [17].

Viscosity of Corn Syrup and Mucus
The viscosity of corn syrup is roughly known to between 2-4 Pa · s, depending on the type of corn syrup. One study provided a value of corn syrup viscosity of ∼3.2 Pa · s, extrapolated for when concentration of corn syrup is 100% [23]. In our study, however, the viscosity of corn syrup was measured to be ∼5.7 Pa · s.
The corn syrup used in this experiment (Kroger Light Corn Syrup-32 fl oz (946 mL), and UPC 011110738639) was stored at room temperature. The viscosity of the corn syrup was determined using a viscosity meter (AnD SV Series Tuning Fork Vibro Viscosity meter, A&D Store, SV-10A). The viscosity meter was calibrated using a two-point calibration method, using distilled water, viscosity 0.0009 Pa · s and glycerol, viscosity 1.137 Pa · s. Temperature of the fluids used in calibration was measured by the viscosity meter and was taken into consideration when determining their known viscosity. Temperature was measured ∼22.4 degrees Celsius during calibration.
Both the samples used for calibration, and the corn syrup, was measured in 10 mL sample containers (AND Instruments AX-SV-34 10 mL Sample Cup), and viscosity was measured using 10 mL of each fluid. Since the viscosity of the corn syrup used in our experiment had a higher than expected value, ten separate trials were done, with different samples, from the same bottle of the Kroger Light Corn Syrup, to verify the viscosity. The resulting average of the ten trials was 5.744 Pa · s, with a standard deviation of 0.7017. The temperature of the corn syrup was also determined by the viscosity meter. The average temperature measured across ten trials was 22.44 degrees Celsius, with a standard deviation of 0.2633. We can conclude the viscosity of our corn syrup is therefore ∼5.7 Pa · s, and is indeed higher than theoretically assumed from prior investigations.
NHTE mucus was a fluid with unknown viscosity that needed to be determined in situ. Calculations of R allow η, the viscosity of mucus, to be estimated. Using Equations (4)-(12), we could express R as a function of η. Then, by plugging in different values for η and calculating the corresponding R, an appropriate value for η could be found. Table 2 shows R for various values of η. We observed, from TEM images in Ref. [17], that R for PVP-FeNps in mucus were around 0.5 to 1 µm. From Table 2, we see that a η of 15 Pa · s would produce the appropriate R. Therefore, we estimate η of mucus to be about 15 Pa · s. Table 2. Viscosity estimates for corresponding MNP aggregate radii in human tracheobronchial epithelial (NHTE) cell mucus when exposed to the magnetic field. For each given value of η, aggregate radii R for polyvinylpyrrolidone coated Fe 3 O 4 nanoparticles (PVP-FeNPs) and poly(lactic-co-glycolic) acid-dexamethasone FITC BaFe 12 O 19 nanoparticles (PLGA-Dex-FITC-BaNPs) were estimated using Equations (4)- (12). This estimate included calculations using r m and r n , with their values and uncertainties listed in Table 1. Uncertainty in R was then propagated from errors in r m and r n by the method of Section 2.8.

Fokker-Planck Simulation and Survival Probability
Magnetic nanoparticle movements were simulated with the one-dimensional, Fokker-Planck equation [24]: where P(z, t) is the one-dimensional probability density function, D is the diffusion coefficient given by the Stokes-Einstein equation [25], Here, k b is the Boltzmann's constant, T is the temperature (300 • K). v(z, t) is the NP drift speed due to external forces. In this case, since all external forces are time-independent and v ter is reached almost instantaneously, the drift speed v(z, t) can be replaced by v ter (z) as a function of position z only. All Fokker-Planck simulations were performed with Mathematica (version 12.1, Wolfram, Champaign, IL, USA), a symbolic mathematical computation program, by solving the Fokker-Planck differential given in Equation (13) for P(z, t) numerically, with appropriate boundary conditions and methods. Specifically, the finite element method of lines was used [26]. For a list of appropriate differential solver methods, refer to Wolfram (Champaign, IL, USA) documentation [27]. Once P(z, t) is found, the survival probability, Q(t), which represents the percentage of MNPs remaining in the fluid over time, can be calculated by where a and b are positions that border the fluid of interest.

Magnetic Field Simulation
The magnetic pulling force produced by the magnetic tip shown in Figure 1 was simulated by equation: where B(z) was obtained from a curve fit shown in Figure 5 with a 9th-degree polynomial, M sat is saturation magnetization, and V is the volume of the NP. Magnetic field setup, data, and saturation magnetization are taken from Ref. [17].

Error Propagation
All error propagations, except for the graphical estimate of v ter , were based on error propagation theory and formulas from Ref. [28]. The propagation of uncertainty is given by where x, ..., z are quantities measured with uncertainties δx, ..., δz, and q is a function of x, ..., z. The formula calculates δq, the expected uncertainty in q. Error propagation was performed in Mathematica.

Results
Graphical results of the Fokker-Planck simulation are presented in this section. First, τ of a typical observed MNP aggregate was calculated. Then, we estimated MNPs v ter using both data from Ref. [17] and the graphical approach described in Section 2.3. These velocities allowed us to calculate aggregate radius R for each MNP. Finally, using calculated M sat and R, we simulated NP movements with Fokker-Planck, Equation (13). After Fokker-Planck Simulations, probability density functions of NPs were presented graphically, and survival probabilities were plotted as a function of time.

Characteristic Time
A typical NP will reach approximately 95% of its v ter within three τs [11]. For MNPs, similar calculations show that 99% of their v ter is reached within a fraction of one τ. Characteristic time τ is a function of η, r m , n m , R, and m, using Equations (2), (6), (11), and (12): Thus, small MNPs in viscous fluids will reach their v ter almost immediately (see Figure 6) because τ is small. Direct observation of v ter is only possible because of the negligible τ. Tables 3 and 4 show values of τ for some select MNPs.

Typical MNP in Gravity
With  (1) and (3), respectively. Note that the displayed velocity for poly(lactic-co-glycolic) acid-dexamethasone FITC BaFe 12 O 19 nanoparticle (PLGA-Dex-FITC-BaNp) in MF is only valid at a position 2500 µm above the magnetic tip, and will increase if the BaNp is positioned closer to the magnetic tip. Table 3. Characteristic time, τ of FeNps in corn syrup (CS) and mucus. All τs are 10 −7 s or less and hence have little effect on v ter . Characteristic time τ was calculated by Equations (7), (9), (10), and (18). Uncertainty in R was propagated from errors in r m and r n listed in Table 1, by the method of Section 2.8. The propagated uncertainties in τ were negligible.

Saturation Magnetization
Milling and Surface modifications are known to decrease saturation magnetization of MNPs. This is likely due to the procedures' destructive effects on the crystal structure of MNPs. Milling breaks big magnetic particles into smaller sizes usable for mucus penetration. With milling, the saturation magnetization of BaNPs drops 18 fold, drastically reducing its magnetic pulling power [17]. In surface modification where the entire surface of MNPs is covered with functional groups, saturation magnetization is reduced by 14.5% [29]. In the study by [17], each surface modification covers the MNP surface by about 10%, thereby reducing its magnetization by approximately 1.45% each time a specific molecule is functionalized to the surface. The results of these estimates are summarized in Table 5. Table 5. Modification/milling effects on saturation magnetization of different MNP types. Each successive surface modification decreases the MNPs saturation magnetization by 1.45%, while milling decreases it by 95%, i.e., DL-2-amino-3-phosphonopropionic acid (AP3) BaNp. Saturation magnetization M sat and its uncertainty was estimated from magnetometry measurements given in Ref. [17] and propagated from errors by the method of Section 2.8 for each surface modification.

Terminal Velocity
In general, v ter increased greatly when a static MF is applied. However, the non-MF particle velocities and the amount of v ter increase in MF also depended on MNP type, surface modification, and the liquid medium. FeNP speeds increased significantly when an MF is applied, while BaNPs had a faster non-MF speed but slower MF speed relative to FeNPs. Surface modifications decreased MNPs speeds, and, in mucus, all MNP speeds are drastically reduced, see Table 6. Table 6. Terminal velocity, v ter of different MNP types and surface modifications, in mucus or corn syrup (CS), with gravity only or with gravity and MF. Terminal velocities are constant under gravity, whereas under a spatial varying MF, v ter s depended on the position z above the magnet. For corn syrup, v ter was recorded at 5000 µm above the magnetic tip, while, for mucus, it was recorded at 2450 µm above the magnetic tip. Terminal velocity v ter was calculated from r m , r n , and M sat with values and errors listed in Tables 1 and 5, by Equations (8), (11), and (12). Uncertainty in v ter was then propagated from errors in r m , r n , and M sat by the method of Section 2.8. We note that v ter could also be determined by the graphical approach in Section 2.3.

Surface Modification
Liquid

MNP Aggregate Radius
MNPs aggregate radii, R m , calculated by Equation (9), are listed in Table 7 for gravity and MF. All surface modifications reduced R m , and it decreased with an increasing number of modifications. All R m s were further reduced in mucus. Table 1 lists single MNP and NMP sizes before they formed aggregates. Table 7. Surface modification effects on R m , with gravity or MF, in corn syrup (CS) and mucus. Magnetic Radius R m was estimated from r m with values and errors listed in Table 1, by Equation (9). Uncertainty in R m was propagated from errors by the method of Section 2.8.

Fokker-Planck Simulation
In the following simulations, MNPs diffuse both horizontally along the x-axis and vertically along the z-axis, with initial probability densities described by Equations (19) and (20), respectively.
These equations describe normalized Dirac-Delta like functions centered at x = 0 µm (average initial horizontal position of MNPs) and z = 7500 µm (average initial vertical position of MNPs; distance above the magnetic tip). Figures 7 and 8 show these initial positions. MNPs are pulled by gravity and a static magnetic field in the negative z-direction.
In the x-direction, no forces act on the particles and hence there is no velocity term, and the Fokker-Planck equation given in Equation (13) reduces to a simple diffusion equation. In viscous environments such as corn syrup (η ∼ 5.7 Pa · s) and NHTE mucus (η ∼ 15 Pa · s), all MNP simulations showed no change in probability density function P(x, t) in the xdirection over 24 h, see Figure 7. In the z-direction where gravity and MF are pulling on MNPs, 24-h simulations also showed no change in P(z, t) when MNPs do not combine to form large aggregates, see Magnetic nanoparticles were only pulled when an MF was acting on their aggregates. Figure 9 shows a simulation for FeNp aggregates under a high gradient MF. Initially, at x = 7500 µm, FeNPs were pulled through a corn syrup layer of 5000 µm thickness to their final positions around x = 2500 µm within a few seconds. Figure 10 shows this simulation in mucus. In mucus, MNPs took approximately 30 min to one hour to cross a 100 µm thick mucus layer. The widths of the probability functions were also much narrower in mucus than in corn syrup. These spikes disappeared when a Mathematica numerical cell size smaller than 1 µm was used but that took longer to simulate.

Survival Probability
Fokker-Planck simulations also revealed that the magnetic force was the major force responsible for the movement of MNPs. Figure 11 shows the survival probability Q(t) calculated by Equation (15) in corn syrup. The survival probability of any MNPs dropped from 100% to zero only when an MF was applied. Figure 12 shows the corresponding survival probability in mucus.  Figure 11. FeNP aggregates' survival probability Q(t) in corn syrup. Blue circles represent particle motion under gravity and diffusion. Orange squares represent particles that are also exposed to a magnetic field.

Discussion
In this study, we were able to simulate MNP movements with the Fokker-Planck equation and estimate MNP aggregate sizes. We found diffusion effects estimated from Stokes' law to be small in mucus, giving a very narrow P(z, t) width. We also found R m to decrease with the number of surface modifications and to decrease with fluids of increasing η.
Deviation from Stokes-Einstein relation was detected when we compared Fokker-Planck simulation in the x-direction (diffusion only) in Figure 7 with MNP behavior observed from Ref. [17] in Figures 13 and 14. Figure 7 showed no change from initial positions over time, but there is irregular spreading in Figures 13 and 14. Furthermore, comparing Figure 13 with electron microscope images of the mucus structure depicted in [30], we see that the PVP-FeNPs are conforming to a particular mucus structure scaffold, rather than undergoing simple Brownian diffusion. The mucus scaffold shows mucus molecules (mucin) to be linked and to form a scaffold structure of unknown rigidity. The rigidity of this structure could hinder the normal diffusion process described by Equation (14). In addition, instead of going through simple Brownian diffusion, MNP particles would be pulled through the mucus pores by the MF. While the particles are being pulled along the negative z-direction, P(x, t) in the x-coordinate could also be affected by that pulling force. Therefore, it is possible that Stokes' law has to be modified, and a 3D simulation needs to be done in this case.  There is no magnetic force and gravity in the horizontal direction, only diffusion effects exist. z = 0 represents the bottom of mucus layer. Images were taken at 10, 20, and 30 min after releasing MNPs. This figure is taken directly from Ref. [17]. Reproduced with permission from Economou, E.C, BioNanoScience; published by Springer Nature, 2016.
Deviation from the Fokker-Planck equation was also detected when we compared Fokker-Planck mucus simulation in the z-direction (diffusion+gravity+MF) in Figures 4 and 10 from Ref. [17]. Figure 4 shows a much wider MNP distribution than those of the mucus simulations. In Figure 4, the initial density function peak is missing at time t = 10 min, and the density is smoothed out all over the mucus layer. Then, near the bottom of the mucus at t = 30 min, MNPs slowed down to form a peak around x = 2 µm rather than settling down to the bottom of the mucus (the expected result in a high MF gradient would have shown little or no peak near the bottom). Again, this deviation is best explained by mucus structure scaffold and aggregate sizes. Near the top of mucus, MNPs aggregate sizes are small relative to pore size but MNP distribution is restricted to a small area in the x-direction, so density function peak was limited by the number of MNPs that can pass through a position z µm above the magnetic tip. Then, near the bottom of the mucus, a high gradient MF causes MNPs to increase in aggregate size. The increased sizes hinder their passage through scaffold pores, but MNP distribution is spread out in the x-direction.
Thus, many MNPs could pass through a position z, but, at a slower speed, as some MNPs are momentarily trapped by the scaffold pores.
Despite this scaffolding structure, large and small MNPs alike were pulled through the mucus layer by MF with surprising efficiency. This might be due in part to the viscoelastic properties of human lung mucus [31]. Under deformation by a magnetic force, the lung pore structures can exhibit both viscous and elastic behavior, allowing MNP clusters with sizes much greater than those of mucus pores to be pulled through. MNP movements through mucus allow determination of η from Table 2, except in the case of milled particles such as PLGA-Dex-FITC-BaNp. At a η of 15 Pa · s, Table 2 shows that PLGA-Dex-FITC-BaNp has a R of ∼ 32 µm. However, TEM images in Ref. [17] show that clusters were around 1µm, corresponding to a viscosity of η = 150 Pa·s according to Table 2. Smaller clusters mean that magnetic force per unit mass was higher, so this discrepancy could have arisen from the estimation of M sat for milled BaNps. We suspect that, even though M sat of BaNps suffered a 95% loss from milling (see Table 5), during the formation of aggregates, some BaNp lattice structures were able to realign themselves, thus recovering some of their magnetic abilities. Therefore, our measurement suggests that NHTE cell mucus viscosity ranges between 15 Pa·s to 150 Pa · s. This result is well within the range of reported human mucus viscosity values that can be as high as 10 4 − 10 6 Pa · s and as low as 10 −2 Pa · s [32].
The ability to estimate both MNP radius and fluid η in situ would prove useful in applications such as alternating current (AC) magnetic hyperthermia [33]. Heating and other physical dynamics of MNPs change according to the size of the aggregate, and MNP aggregates need to be within a specific size range for some applications to be effective. In addition, viscosity η, a measure of a substance's resistance to the motion under an applied force, is generally unknown in various biological fluids. Using methods like that of Table 2, the range for η could be estimated. Being able to find the η range of a biological fluid will allow better design of MNPs in future applications, such as estimating the required magnetic field strength to pull MNP through high viscosity fluids [34].
However, to improve the accuracy of these estimates, Stokes' law, Stokes-Einstein relation, or Fokker-Planck equation must be modified. The drag coefficient from Stokes' law (Equations (5) and (6)) appears in the Fokker-Planck equation (13) twice. First, in the Stokes-Einstein relation, Equation (14), then, in the terminal velocity term, see Equation (8). We saw that, in the x-direction, particle diffusion in experiments does not match Stokes-Einstein's relation, and the MNP distribution appears different for different MNPs, perhaps due to their specific interactions with the mucus scaffold. We also saw that, in the zdirection, particle distribution overtime does not match P(z, t) simulations. Therefore, v ter from Equation (8) can only serve as a rough estimate unless the drag force in Stokes' law is modified, Stokes-Einstein's relation is replaced, or the Fokker-Planck equation is somehow corrected.
To model the mucus scaffold in future simulations, the drag force could include two terms, size drag, and density drag, defined by and where F size and F density are size and density drag forces, respectively, R m is the radius of MNP aggregates, and R p is the radius of mucus pores. An increase in the MNP aggregate radius would increase the drag force, while an increase in pore radius decreases the drag force. Here, R m is also a function of the MF strength, since stronger MF gives larger R m . Similarly, ρ s and ρ m are densities of mucus scaffold and MNP aggregates, respectively, and the density drag force is defined by their ratio. A and B are proportionality constants that differ between different MNPs because each of them has different chemical interactions with the scaffold. In addition, to account for the number of MNPs that can pass through mucus pores at various distances z from the tip of the magnet, P(x, t) along the x-axis must be considered when simulating P(z, t) along the z-axis. Therefore, it might be best to simulate the Fokker-Planck equation in 3D, rather than x and z in 1D separately. Additionally, some chemical aspects made it hard to estimate sizes and velocities accurately with this model. According to experiments performed by Economou et. al. [17], many iron oxide nanoparticles are seen attached to a huge PLGA-Dex globule (NMP) when the MNPs were generated. If the number of MNPs and NMPs were to increase, then it seems reasonable to have a 100 MNP to 1 NMP ratio in an aggregate. However, when PLGA concentrations were reduced during the generation process [17], the resulting PLGA-Dex-BaNPs v ter greatly increased. Reducing PLGA concentrations likely resulted in a decreased PLGA number in an aggregate. Reducing PLGA number per aggregate would drastically decrease the NP size, causing drag force to decrease and the aggregate v ter to increase. Therefore, the exact MNP to NMP ratio would depend on the concentration of NMPs used when creating these NPs, and cannot be easily estimated without first knowing how solute concentrations affect the MNP to NMP ratio during the NP generation process.

Conclusions
We were able to discover some properties of magnetic nanoparticles from our research. First, there was a clear relation between surface modification and magnetic nanoparticle's aggregate size. Then, we found that iron oxide, Fe 3 O 4 , nanoparticles (FeNPs) and barium hexaferrite, BaFe 12 O 19 , and nanoparticles (BaNPs) react differently to an applied magnetic field, and from iron oxide nanoparticle aggregate sizes, we could estimate the viscosity of the surrounding fluid. However, we also discovered that Stokes' law, Stokes-Einstein's relation, or Fokker-Planck simulations do not match magnetic nanoparticle behavior in the selected experiments performed in a mucus layer of air-liquid interface cultured primary normal human tracheobronchial epithelial cells by Economou et al. [17]. Thus, future modifications to these equations must be considered. In particular, modeling the mucus scaffold might require a very different drag force than that of spheres falling through laminar fluid flow in Stokes' law.