Imaging of pH distribution inside individual microdroplet by stimulated Raman microscopy

As the dominant inorganic pair of an acid–conjugate base in ambient aerosols, sulfate (SO₄²⁻) and bisulfate (HSO₄⁻) have distinguishable vibrational shifts for ν_s(SO₄²⁻) and ν_s(HSO₄⁻) at ~985 and ~1,050 cm⁻¹, respectively (Fig. 1B). We first measured the SRS spectra of Na₂SO₄ and NaHSO₄ in standard bulk buffer solutions to establish calibration curves of [SO₄²⁻] and [HSO₄⁻] (SI Appendix, Fig. S1). As expected, the intensity of ν_s(SO₄²⁻) and ν_s(HSO₄⁻) increased with concentration (Fig. 1B), with a good linear correlation (SI Appendix, Fig. S2). Meanwhile, the intensity of ν(H₂O) (3,420 cm⁻¹) maintains relatively stable (Fig. 1C). Geometric factors of microdroplets and the wave nature of light are the major factors that hinder the data interpretation (25); hence, the internal standards (IS) are adapted to correct these possible influences (SI Appendix, Note S2). The intensity ratio of ν_s(SO₄²⁻)/v(H₂O) and ν_s(HSO₄⁻)/v(H₂O) exhibits an excellent linear relationship with the [SO₄²⁻] and [HSO₄⁻] (R² = 1.000 and 0.994, Fig. 1D), which allows the quantification of [SO₄²⁻] and [HSO₄⁻]. We also showcased the reliability of this linear fit by repeating this IS method with the MRS (13, 25, 26) (SI Appendix, Figs. S3 and S4). Although the detected [SO₄²⁻]_MRS and [HSO₄⁻]_MRS through MRS are similar to those of SRS (Fig. 1E), the latter exhibits significantly higher precision due to a stable calibration curve (detailed in SI Appendix, Note S3). The sum of the detected [SO₄²⁻]_SRS and [HSO₄⁻]_SRS of a series of mixed solutions is almost equal to the preknown total [S] with a P-value of 0.941. These consistencies solidify the potential of SRS for pH detection by high discrimination for the dissociation of SO₄²⁻/HSO₄⁻.

To elucidate the pH distribution within whole microdroplet aerosol, we have derived an analytical expression for pH calculation with [SO₄²⁻] and [HSO₄⁻] based on the definition of acidity and Debye–Hückel theory (27), as detailed in SI Appendix, Note S4.

The pK_a is the dissociation constant of HSO₄⁻, considered as a fixed value (a reasonable explanation is provided in SI Appendix, Note S5). Log(γ_bisulfate/γ_sulfate), a derived term, may represent the effect of ionic strength on the activity coefficient. Physically, the activity coefficient of ions comes from the effect of interionic forces and ion-solvent molecular forces on the free energy of ions in solution (28). As shown in Fig. 1F, log(γ_bisulfate/γ_sulfate) exhibits an obvious linear relationship with $\sqrt{4\times [{\mathrm{SO}}_4^{2-}]+[{\mathrm{HSO}}_4^{-}]}$ that derived from the ionic strength formula. Then, log(γ_bisulfate/γ_sulfate) can be obtained via Eq. 1B.

where a and b, the constants of the calibration curve, were determined at 0.01164 and 0.4748 by the linear fitting (Fig. 1F), respectively, wherein the log(γ_bisulfate/γ_sulfate) was solved with Eq. 1A and preknown pH of NaHSO₄ standard solutions and mixed solutions. Despite slight fluctuations, the difference between the estimated pH (via Eqs. 1A and 1B) and preknown pH in Fig. 1F is insignificant (P = 0.976, one-tailed t test). Meanwhile, the P-value of estimated pH via MRS and preknown pH is 0.054 (SI Appendix, Fig. S5), indicating the accuracy of SRS and versatility of Eqs. 1A and 1B. Notably, the regressive Eq. 1B cannot cover the high-concentration cases within the supersaturated aerosol, but the thermodynamic model evaluation demonstrated the practicability of Eq. 1B for estimating the pH of microdroplet with high log(γ_bisulfate/γ_sulfate) (SI Appendix, Note S5 and Figs. S6 and S7).

Aerosol microdroplet with a diameter of 2.9 μm was characterized to explore the possible two-dimensional (2D) pH distribution within their bulk under a steady-state environment based on a superhydrophobic substrate (SI Appendix, Note S6 and Figs. S8 and S9 and Movie S1). Fig. 2A shows SRS images that track the Raman intensities of ν(H₂O), ν_s(SO₄²⁻), and ν_s(HSO₄⁻) of a ~2.9-μm microdroplet. The SRS image of ν(H₂O) demonstrated a higher intensity in the centroid region, which decreased as the objective moved to the microdroplet edge owing to the reduced actual volume that was probed by SRS. Unexpectedly, the intensities of both ν_s(SO₄²⁻) and ν_s(HSO₄⁻) were higher at the microdroplet edge. Fig. 2B plotted the distributions of [SO₄²⁻], [HSO₄⁻], as well as pH within a microdroplet. Evidently, all these substances including [H⁺] (10^−pH) exhibited a strong tendency of interfacial enrichment, with a sharp increase of concentrations from the microdroplet center to the edge. This concentration gradient within a microdroplet was confirmed by MRS as well (SI Appendix, Fig. S10), despite the lower spatial resolution. Similarly, for bulk solutions, the species at the air–water interface also represented enrichment (SI Appendix, Fig. S11), consistent with the trend observed in microdroplets. The concentration gradient of [SO₄²⁻] diminished as the probe laser moved deep into the bulk phase of the solution (SI Appendix, Fig. S12).

The distributions of [SO₄²⁻], [HSO₄⁻], and pH within a 2.9-μm microdroplet were detailed via the line profile and its representative average (SI Appendix, Note S7 and Fig. S13). As shown in Fig. 2C, the distributions of [SO₄²⁻] and [HSO₄⁻] derived from four line profiles (marked by purple dashed lines in the inset) exhibited very similar trends, implying that these distributions are centrosymmetric. Specifically, the [SO₄²⁻] increased from 1803.5 to 3306.2 mM (1.8-fold increase) as the probe laser moved 1.3-μm away from the microdroplet center, and then cannot be detected when the probe moved further away. Meanwhile, [HSO₄⁻] increased from 2746.6 to 9443.8 mM (3.4-fold increase) from 0.0 (center) to 1.5 μm (edge), indicating a stronger interfacial enrichment (Fig. 2 C, Inset). In contrast to the strong enrichment of HSO₄⁻ at the interface, [SO₄²⁻] may be rapidly reduced below the detection limit of SRS due to the rapid conversion of SO₄²⁻ to HSO₄⁻ under low pH (<−0.41). Based on the distributions of [SO₄²⁻] and [HSO₄⁻], we found that the pH within the 2.9-μm microdroplet gradually decreased from center to edge. The pH (−0.34) at the microdroplet edge is much more acidic than the centroid (0.19), with [H⁺] increased to 2.19 M at the edge (3.4-fold compared with the centroid).

This pH distribution was further examined through classical molecular dynamics (MD) simulations (29). Fig. 2D depicts the time-dependent snapshots of H₃O⁺, HSO₄⁻, Na⁺, and SO₄²⁻ in a water nanodroplet (10 nm) from 0 to 50 ns. These species migrate toward the droplet edge, achieving a relatively stable distribution after 40 ns. Fig. 2E shows the total number profiles of these species along the radius direction (from center to edge) at 50 ns. Species number within the nanodroplet obviously increased from center to edge, which is in accordance with the experimental results (Fig. 2 A–C). Besides, the number density distribution of H₃O⁺, HSO₄⁻, and SO₄²⁻ within a nanodroplet (Fig. 2F) is comparable with the chemical imaging results of a small microdroplet (Fig. 2B). Here, simulated H₃O⁺ ions have a relatively greater propensity for the air–water interface compared with other ions, supporting the experimental results of the detected pH distribution in the aerosol microdroplets.

Further, the 3D chemical distribution of a bigger microdroplet (27.0 μm) with a lower sulfates concentration was reconstructed (Movie S2) to explore the pH distribution within atmospheric cloud microdroplets. Fig. 3A shows its 2D pH distribution at different Z values (SO₄²⁻ or HSO₄⁻ distribution is in SI Appendix, Fig. S14). The distribution of pH and (bi)sulfate at different Z-values is central symmetric, with lower pH at edges in any direction compared with the microdroplet center. Meanwhile, Fig. 3B shows that the pH distribution in the XZ plane is axisymmetric and similar to the one in the XY plane at Z = 0.0 μm (white box in Fig. 2D). These results demonstrated a central symmetric pH distribution within the entire volume of cloud microdroplets.

The Fig. 3C details this pH distribution where pH remains stable (~0.27) within the distance of 9.5 µm from the center, and down to a pH value of ~ −0.09 at the edge. The difference (0.69 M) in [H⁺] between the center and edge here is much smaller than that (1.54 M) for a 2.9-μm microdroplet, indicating that the particle size plays a key role in the H⁺ enrichment at the interface. Interestingly, this process of pH decreasing from center to edge is not always monotonous. There is an unexpected pH increase at the inner layer near the edge of the 27.0-µm microdroplet, which is not observed in the case of the 2.9-μm microdroplet. In detail, from the center toward the microdroplet edge, along X, Y, and Z directions, pH first increases from 0.27 to 0.35, then sharply decreases to −0.09 at the interface. Fig. 3D shows the constructed model of 3D pH distribution in small aerosol and big cloud microdroplets based on the aqueous microdroplets with sizes of 2.9 and 27.0 μm, respectively. The notable difference is that cloud microdroplets possess a nonmonotonic pH variation and an abnormal increase in pH (green region) occurs near the edge of larger droplets. Lhee et al. likewise found that the concentration of charged dyes progressively increased from the center to the edge of 7 μm droplets, while a decrease of the concentration was observed close to the interface in the cases of 24 and 28 μm droplets (30). This phenomenon in larger microdroplets can be attributed to the migration of molecules from the adjacent inner region to the interface, which may be driven by the interfacial electric field (EF) (9).

The Fig. 4A details the effect of droplet size on pH distribution as well as acidity change (ΔpH) of microdroplets. The results show pH variation inside typical aerosol microdroplets as a function of normalized distance to the center. Overall, the smaller the size of the droplet is, the lower its interfacial pH and the greater the degree and zone of pH variation are. For small microdroplets (d ≤ 3.2 μm), the pH decreases continuously from the center (0.11 to 0.37) to the edge (−0.35 to −0.78). Whereas, as the diameter increases, a “pH-stable zone” begins to appear in the central region of the microdroplet while the edges remain more acidic. The width of this pH-stable zone increases with increasing droplet size; in other words, the “acidification zone” (the zone with the big absolute value of ΔpH) decreases accordingly. For example, the width of the acidification zone for 134.6 μm microdroplets is 0.137 (normalized distance, ND), which is almost 6 times smaller than that of the 2.9-μm microdroplet (ND = ~0.8). Additionally, for small microdroplets (d ≤ 3.2 μm), they showed a monotonic trend of pH decreasing from center to edge. For bigger microdroplets (21.9 μm ≤ d ≤ 134.6 μm), they all presented a sharp pH decrease in the region close to the interface, while the “pH-stable zone” is not entirely similar. When the size is larger than 21.9 μm, the pH first slowly increased from the center to the inner region of the edge, then sharply decreased to the lowest pH at the edge. The highest pH is observed at the inner region of the edge. Based on the statistical H⁺ distribution of five microdroplets of different sizes, we proposed a potential method to estimate the average pH at the interface of ambient aerosols (SI Appendix, Note S8 and Fig. S15).

As reported, ions are dragged to the interface from the bulk through the electrostatic interactions and accumulate in nearby inner layers (31). Therefore, the [H⁺] at the microdroplet edge could partially reflect the electric field strength. As shown in Fig. 4B, the average edge pH possesses a good linear relationship with Δ[H⁺], demonstrating the possible influence of average [H⁺] at the microdroplet edge on the pH gradient within the microdroplet. Moreover, the distribution of simulated EF (SI Appendix, Fig. S16) in microdroplets is consistent with that of species (32). The EF gradients are much larger at the edge of big microdroplets than that of small microdroplets (Fig. 4C). It implies that this pH gradient may be caused by the electrostatic interaction force from the interfacial electric field. This force migrates H⁺ from the inside of microdroplets to its edge. Then, sulfates are protonated with H⁺ leading to the enrichment of bisulfates. Fig. 4D shows the relationship between interfacial affinities of (bi)sulfate and microdroplet size. Obviously, SO₄²⁻ preferentially resides deeper in the interfacial region relative to HSO₄⁻, and only HSO₄⁻ was detected in the microdroplet edge, consistent with the reported results (33). The normalized width layer containing only HSO₄⁻ (green columns) decreases from 0.12 to 0.02 when the microdroplet diameter increases from 2.5 to 33.9 μm, then remains stable (~0.01) within the bigger microdroplet. The surface preference of HSO₄⁻ becomes more obvious in small droplets, following the abovementioned lower pH and enhanced electric field at the surface of small microdroplets (detailed in SI Appendix, Note S9).

As illustrated in Fig. 4E, the Δ[H⁺] is proportional to S/V with a surprising R² of 0.94. Meanwhile, the width (ΔL) of the acidification zone has a linear relationship with the logarithm of S/V. Accordingly, the pH distribution inside microdroplets can be estimated by S/V, which is favorable for improving the accuracy of the aerosol model prediction. In addition, we also investigated the pH distribution under low relative humidity (RH) of 35%. Compared to microdroplets under 90% RH (Fig. 3A), the pH gradient within the microdroplet under 35% RH was still present but less pronounced (SI Appendix, Fig. S17), especially in the one with a diameter of 9.0 μm. For instance, at 35% RH, the [H⁺]_edge/[H⁺]_center of microdroplets with diameters of 9.0 and 24.4 μm is 2.26 and 1.64, respectively, while under 90% RH with a similar diameter, this ratio can reach 3.78 (9.1 μm-microdroplet) and 2.02 (27.0 μm-microdroplet) (Fig. 4A), respectively. The increased viscosity in the concentrated aerosol may inhibit ion migration and electric field formation (34, 35), resulting in a weaker pH gradient. For big microdroplets with radii of 23.6 and 38.2 μm, they are more prone to crystallize in the center of the microdroplet, leading to phase separation. The residual liquid phase of the microdroplet aerosol still exhibits a pH gradient and possesses a more acidic interface.

Like our previous work (44), the MRS of bulk solutions and aerosols in a quartz cell (SI Appendix, Fig. S1A) was recorded by XploRA Plus confocal spectrometer (Jobin Yvon, Horiba Gr, France) with 100× Olympus microscope objective (LMPLFLN100×, Olympus, N.A. = 0.8). The power of the external-cavity continuous-wave diode laser (532 nm) was set to 45 mW (50% of the maximum power) to improve the signal-to-noise ratio. Before the measurement, the instrument was calibrated against the Stokes Raman signal of pure Si [a silicon wafer ([110] crystal plane)] at 520 cm⁻¹. Every spectral resolution was 1 cm⁻¹ over the range of 300 to 4,500 cm⁻¹ with an acquisition time of 10 s and three accumulations. To obtain the 2D SO₄²⁻/HSO₄⁻ and pH distribution of microdroplets, Raman mapping was conducted based on a motorized scanning table with a minimum step size of 0.2 μm. Each pixel represents a single Raman spectrum collected with the same setup as a bulk solution.

Our chemical imaging of SO₄²⁻, HSO₄⁻, and H₂O was based on a conventional SRS microscope as detailed in our previous work (22, 45). In brief, a commercial femtosecond laser system (Insight DS+, Spectra-Physics) produced two synchronized pulse trains at 80MHz. The fixed fundamental output of 1,040 nm was employed as the Stokes beam, while the tunable optical parametric oscillator output (680 to 1,300 nm) served as the pump beam. To acquire high spectral resolution (~13 cm⁻¹), pulse durations of the pump and Stokes beams were chirped and stretched to several picoseconds by passing through SF57 glass rods. The intensity of the 1,040 nm Stokes beam was modulated at 1/4 of the laser pulse repetition rate using an electro-optical modulator (EOM, EO-AM-R-20-C2, Thorlabs). The two laser beams were spatially and temporally overlapped via a dichroic mirror, then delivered into a laser scanning microscope (FV1200, Olympus) equipped with galvo mirrors for raster scanning. The combined beams were focused onto the microdroplet by a 60× water immersion objective lens (Olympus, UPLSAPO 60XWIR, NA 1.2). The stimulated Raman loss signal was optically filtered (CARS ET890/220, Chroma), detected by a homemade back-biased photodiode (PD), and demodulated with a lock-in amplifier (HF2LI, Zurich Instruments) to feed the analog input of the microscope to form images with the size of 512 × 512 pixels. The microdroplet was kept in a closed homemade chamber to ensure stable chemical imaging. The laser powers of 20 mW for the pump and 30 mW for Stokes pulses were used for samples. A lock-in time constant was set as 2 to 10 μs to match the pixel dwell time. All images were taken in transmission mode. The Raman bands included the ν_s(SO₄²⁻) and the ν_s(HSO₄⁻) which were recorded at the pump laser wavelength of 940 nm, while the signal of v(H₂O) was measured with the tuned 776 nm pump beam. The optical layout is illustrated in SI Appendix, Fig. S1B.

According to the IS method, standard solutions of Na₂SO₄ and NaHSO₄ within the widest possible range of concentrations are used to create calibration curves relating [SO₄²⁻] and [HSO₄⁻] to the SRS intensity ratio of the ν_s(SO₄²⁻)/v (H₂O) and ν_s(HSO₄⁻)/v(H₂O). The calibration curve of [SO₄²⁻] is established first by a series of Na₂SO₄ solutions (Fig. 1B). After that, the sulfate calibration curve is used to calculate [HSO₄⁻] in the SO₄²⁻/HSO₄⁻ equilibrium of NaHSO₄ standard solutions through elemental conservation (Fig. 1B), obtaining the calibration curve of [HSO₄⁻]. To verify the accuracy of SRS quantitative analysis, similar calibration curves relating [SO₄²⁻] and [HSO₄⁻] to the MRS integrated peak area ratio of the ν_s(SO₄²⁻)/v (H₂O) and ν_s(HSO₄⁻)/v(H₂O) were created first (SI Appendix, Fig. S3). Then, the MRS and SRS spectra (SI Appendix, Fig. S4) of a series of mixed solutions (Fig. 1E) with preknown total [S] were employed to compare the accuracy of two calibration curves: integrated MRS peak area ratio versus concentration (SI Appendix, Fig. S3C) and SRS peak intensity ratio versus concentration (Fig. 1D). All chemicals are purchased from Sigma-Aldrich (Shanghai) Trading Co. Ltd. with high purity (>99%), and the solution pH is detected by a probe (S210, Mettler Toledo).

Microdroplets are generated from the mixed solution of 300 mM NaHSO₄ and 50 mM Na₂SO₄ via a homemade atomizer of 100 mL at room temperature and naturally dispersed on superhydrophobic coverslips. Then, that hydrophobic coverslip with microdroplets is quickly flipped over onto the central perforated slides, forming a closed chamber with the help of another original coverslip and a sealant, as shown in SI Appendix, Fig. S8. Before being probed by SRS or MRS, the chamber was held for 5 min at room condition (23 °C, 101.325 kPa) to ensure that the droplets became stable. The RH in the chamber is 90 ± 2.5%. The pictured micrometer-sized droplets indicated a successful microdroplet collection. Owing to the limitation of the analyzer (SRS/MRS), we prepared the superhydrophobic coverslips as reported in literature (46).

For the distribution of pH, [SO₄²⁻], and [HSO₄⁻], all the images of species distribution were plotted based on corresponding SRS images and MATLAB software. In detail, the SRS image, including H₂O (3,420 cm⁻¹), SO₄²⁻ (985 cm⁻¹), and HSO₄⁻ (1,050 cm⁻¹), in Tag Image File Format, is imported into MATLAB software as a 2D matrix. Consequently, using the cells in the matrix as calculation units, the [SO₄²⁻] and [HSO₄⁻] in any cells can be calculated from the corresponding calibration curves based on SRS. Finally, the figure of (bi)sulfates distribution can be obtained by converting their two-dimensional matrices into images again. Furthermore, the matrices of calculated [SO₄²⁻] and [HSO₄⁻] were used to produce the matrix of pH based on Fig. 1D and Eqs. 1A and 1B, as well as the image of pH distribution. For the line profile, all the images of species distribution were plotted based on the corresponding SRS images and the software of Image J. In detail, the SRS images of H₂O, SO₄²⁻, and HSO₄⁻ were merged into one image by Image J software. Then, the Image J was used to read data of the SRS image in a certain direction. As a result, we obtained line profiles of H₂O, SO₄²⁻, and HSO₄⁻ in one direction, that derived from the same position of the microdroplet (SI Appendix, Fig. S13). After that, the [SO₄²⁻], [HSO₄⁻], and pH can be determined through calibration curves in Fig. 1D and Eqs. 1A and 1B.

To obtain thermodynamic properties and species distribution within microdroplets, a classical MD calculation was performed for sulfate microdroplet systems. An aerosol microdroplet of Na₂SO₄ and NaHSO₄ was constructed in this work. The radius of the aerosol microdroplet is 5 nm, and a 12 × 12 × 12 nm³ simulation box which contains H₃O⁺, Na⁺, HSO₄⁻, SO₄²⁻, and H₂O molecules was used. According to our experimental results, [HSO₄⁻] and [SO₄²⁻] were set at 3.8 and 2.3 M, respectively. The concentration of free H₃O⁺ (activity of the H⁺ ion) was determined by our empirical formula (Eqs. 1A and 1B) to be 1.047 M. Then, the [Na⁺] was determined at 7.353 M based on the conservation of charge. Accordingly, the [H₂O] was calculated at ~36.44 M (29). Periodic boundary conditions were applied in all the three directions (47). The ions were parameterized with the program acpype.py using Antechamber and the GAFF forcefield (48–50). This GAFF forcefield has been widely used for small molecules in atmospheric processes (51, 52). Water molecules were simulated using the transferable intermolecular potential with the three points (TIP3P) model (50). A cutoff distance of 12.0 Å was used for Lennard–Jones and real space coulombic interactions. The Particle Mesh Ewald method was also employed with an interpolation order of 6 with 1.0 Å grid spacing (53). Before the MD began, energy minimizations were carried out. Then, the systems were equilibrated for 200 ps at the temperature of 298 K (NVT). After equilibration, the simulated systems were run for 50 ns in an NVT ensemble at 298 K. The timestep was 1 fs, and the trajectories were collected every 1 ps. The MD trajectories resulting from the dynamic simulations were processed to extract the structural and dynamical properties of the simulated systems. To identify the distribution of ions, the number of ions per 0.5 Å from the center to the air–water interface and the X–Y planar density maps (cross-section of the aerosol microdroplet between 5.5 and 6.5 nm along z-coordinate) was calculated. All the MD calculations were executed by the software package of GROMACS 2019.6. And MD configurations were visualized by the VMD package (54).

M.J., L.Z., and J.S.F. designed research; K.G., and J.A. performed research; L.L. K.G., contributed new reagents/analytic tools; K.G., J.A., K.L., Y.L., G.X., T.W., Z.W., X.Z., L.Z., and J.S.F. analyzed data; K.L., L.L., Y.L., G.X., T.W., and H.C. visualization; H.W. and C.G. review & editing; A.M., H.H., L.W., and J.C. review and editing; and K.G., J.A., H.C., H.W., C.G., A.M., H.H., L.W., L.Z., J.C., and M.J. wrote the paper.

The authors declare no competing interest.