Mechanism of silica–lysozyme composite formation unravelled by in situ fast SAXS

Introduction

A mechanistic understanding of aggregation in aqueous media leading to the formation of composites of inorganic nanoparticles and proteins is of paramount interest for colloid chemistry, Earth sciences, or the design of protein-sensitized biomedical devices and sensors [1-6]. In a broader perspective, many processes involving the assembly of nanoparticles to higher-level hierarchical structures are known to be directed by small organic molecules and macromolecules [7-10]. Such fundamental interactions are the key to understand the crystallisation of biominerals in living organisms (e.g., in bone formation), and to manufacture better functional materials [11-15].

In particular, composites of amorphous silica (SiO₂) nanoparticles (NPs) and lysozyme (LZM) have attracted attention because silica NPs readily form in many Earth surface environments (e.g., oceans, hot springs) where biological activity dominates (e.g., diatom formation) but they are also key components in numerous technological applications from electronics to paint production. In turn, lysozyme adsorbs strongly to silica NPs [16-18], while at the same time preserving its enzymatic activity, most notably antibacterial properties [19]. Over a wide range of pH values (2 to ca. 10), the surface of silica NPs is negatively charged, whereas the LZM molecule is positively charged [20,21]. This way LZM can act as a bridge between silica NPs, leading to aggregation and flocculation and thus to large silica NP–LZM composites. In the SiO₂–LZM model system, a number of studies investigating the relationship between silica NP sizes, and adsorption modes of lysozyme revealed a correlation between composite properties and the folding/structure of the protein, its enzymatic activity and the exact protein localization with respect to the silica NPs inside flocculated composites [16,17,19,21-24]. The inherent nature of the aggregation processes leading to silica–lysozyme composites, involves structural changes at length scales from a few to hundreds of nanometres, which makes this system well-suited to be characterized by scattering methods. In particular, recent advances [25,26] in detector technology for synchrotron-based in situ and time-resolved small-angle X-ray scattering (SAXS), now allow one to follow all steps in the formation of SiO₂–LZM composites from the individual components. Nevertheless, although such scattering data will reflect the in situ state of a system during measurement, the quantitative information related to any changes in structural properties of the particles/species of interest can only be accessed by developing, testing and validating relevant models and bespoke data analysis methods.

In this study we show step-by-step how a scattering model was developed, verified and applied to time-resolved synchrotron-based SAXS data in which we followed in situ the lysozyme-induced aggregation of silica NPs (ca. 5 nm in diameter) at a time resolution of 50 ms. It is important to note that at typical NPs sizes and concentrations [17,24] and protein concentrations [27] usually employed to make such SiO₂–LZM composites, the relative X-ray scattering contrast/intensity of lysozyme, in the presence of silica, is insufficient to determine the protein contribution to the overall scattering pattern (ratio of ca. 1:100 between LZM and NPs), and only the NP component can directly be followed. However, through this study we demonstrate that if the size of the silica NPs is matched with the dimensions of lysozyme (ellipsoidal molecule [28] 3 × 3 × 4.5 nm³), the evolving scattering patterns contain a unique structure-factor contribution originating from the presence of lysozyme, and this way this important contribution can be assessed. The analysis of this structure function through the derived model then allowed us to extract detailed structural information on the deformation of the LZM molecules upon aggregation, and to determine the mechanisms of the formation of SiO₂–LZM composites.

Results and Discussion

This section of the manuscript is organised in the following way: (1) data presentation and a discussion of general trends; (2) derivation of a self-consistent mathematical model to fit the timer-resolved scattering curves; (3) presentation of the quantitative results from the model fits; (4) discussion of the quantitative trends.

Evolution of SAXS patterns and derived aggregation stages

Upon mixing of the silica NPs and the LZM solution we observed very fast flocculation, which indicated the formation of the composites. In Figure 1 we show an overview of these formation processes based on SAXS data collected at a time resolution of 50 ms and spanning ca. 300 s.

In a contour plot of the time-resolved scattering patterns (Figure 1A), one can distinguish four characteristic time periods (I–IV) and one region of interest (ROI V), which spanned through periods II to IV. Period I corresponds to the initial ca. 24 s of the scattering patterns of silica NPs before the injection of lysozyme. Based on this data we determined the initial form factor (size distribution) of the silica NPs prior to mixing with lysozyme (Figure 1B). The scattering pattern in a log–log representation prominently flattens out at low q (i.e., ). This shows that the initial silica NPs were not aggregated and well-suspended. We derived a discrete size distribution (histogram in the inset of Figure 1B) for the NPs from the Monte Carlo fit implemented [29,30] in MCSAS under the a priori assumption that the NPs were spherical in shape [31] (physicochemical parameters of amorphous silica given in Table S1, Supporting Information File 1). The as-obtained histogram indicated that the size distribution was relatively narrow with a mean radius of 2.53 ± 0.01 nm (distribution statistics given in Table S2, Supporting Information File 1). The total integrated volume fraction for the NPs obtained from the fit, was 0.040% ± 0.001%, which matches very closely the expected value of 0.041% calculated for silica NP precipitated from a 1000 ppm SiO₂ solution at pH 7.5 and 21 °C (Table S1, Supporting Information File 1, calculated with PHREEQC [32]).

Period II (ca. 25–30 s) in Figure 1A represents scattering patterns during and soon after the injection of the LZM solution and its mixing with the silica NPs. Period II is hence preceded by a dead-time period of 500 ms (see Experimental section). Period II (Figure 1C) can be divided into multiple steps. The first 1.25 s were primarily characterized by a rapid and significant (ca. 11-fold) increase in intensity at low q (q < 0.3 nm⁻¹). During the following 2 s, the low-q part still kept increasing (to ca. 15-times of the initial intensity) but less rapidly, and at q ≈ 1 nm⁻¹ a characteristic local maximum developed (ROI V in Figure 1A). The intensity increase at low q originated from the formation of large aggregates constituting the composites, with sizes outside the minimum q-range, whereas the local maximum (the correlation peak q ≈ 1 nm⁻¹) indicated the presence of interparticle correlations within those aggregates.

The intensity increase at low q associated with the aggregation continued throughout period III (between ca. 30 and 150 s), yet the correlation peak in ROI V did not change significantly (Figure 1D). Note that up to 150 s (periods I–III) the high-q part of the data (q > 1.5 nm⁻¹, Figure 1C,D) did not change, indicating that the original form factor of silica NPs remained the same after the injection of lysozyme. Therefore, as a first approximation the observed electron-density scattering contrast in these periods (I–III) can be interpreted to originate solely from the silica NPs and not from the combination of silica and lysozyme (Figure S1, Supporting Information File 1). Hence, we could treat the system as a system with two different electron densities (silica NPs and solvent matrix). However, because our silica nanoparticles were smaller (ca. 5 nm) compared to silica NPs in previous scattering studies (ca. 20 nm) on silica–protein composite formation [16,17,23,24] and because the lysozyme addition dramatically changed the silica aggregation state, the contributions of the lysozyme scattering can be accounted for indirectly from the interparticle correlations observed in ROI V. This is key here, because it allows us to extract the changes in lysozyme structural properties as the composites evolve over time, without deriving circumstantial models for a system with three electron densities (silica, lysozyme, and the solvent matrix).

In period IV (150–300 s, Figure 1A and Figure 1D), we observed a further intensity increase at low q (three-times higher at 300 s than at 150 s), which indicated a continuous increase in aggregate size from periods II and III. In time period IV the intensity of ROI V (Figure 1A) started to increase together with the silica form factor at high q (q > 1.5 nm⁻¹, Figure 1D). This suggests that as aggregation continued between 150 and 300 s, the silica NPs themselves started to grow, e.g., due the coalescence of NPs or similar processes (under an assumption that the particles remained spherical in shape).

Scattering model

In the considered silica–LZM composites, the scattering contrast originating from the lysozyme itself can be mostly disregarded (see Figure S1, Supporting Information File 1). Hence, the scattering contrast of the primary silica particles, (Δρ_p)², – expressed as the difference in scattering length density, SLD – is equal to the squared difference in the SLDs of silica and the surrounding water matrix (Table S1, Supporting Information File 1). Consequently, using such an approximation allows us to deal with a two-electron-density system, where the scattering intensity, I(q) is a product of the scattering contrast, (Δρ_p)², the form factor of the silica NPs, P_p(q), weighted by a volume fraction of silica NPs, , and their volume, V_p, together with an effective structure-factor function describing the spatial arrangement of silica NPs within the aggregates, S_eff(q) (Equation 1). We use subscript “p” to emphasize that the parameters and functions concern the primary silica NPs.

Equation 1 is valid only for a system of ideally monodisperse particles (i.e., the distribution is a delta function), which is actually not the case for the silica NPs used here (Figure 1B). This is an important consideration when including the interparticle interactions from the structure factor. There are several approaches to consider the polydispersity of particles together with a structure factor [33], but because the fitted size distribution (histogram in Figure 1B) is discrete with a finite number of n bins a local monodisperse approximation (LMA) [34] is used in our models:

where P_p(q,r_i) is the form factor of a sphere of radius r_i.

The discrete size distribution (Figure 1B) has the form of Equation 3, where for each size contribution (r_i) the corresponding volume fractions, , are known:

Since the partial i-th structure factor, S_eff(q,r_i), includes the interparticle correlations between silica NPs and lysozyme (the local maximum at q_max ≈ 1 nm⁻¹, ROI V in Figure 1), it is also dependent on r_i of the primary silica NPs. Additionally, the structure factor expression has to account for the aggregation of the silica NPs to large objects (the low-q increase), yet the size of these aggregates, in turn, does not necessarily depend on the size of the primary silica NPs.

The interparticle correlations and the local maximum

Under the considered physicochemical conditions the inorganic silica NPs and the protein molecules are oppositely charged [21], and hence they interact through the attractive potential. This, in turn, leads to the formation of aggregates/composites in which NPs behave closely to adhesive hard spheres. The aforementioned interactions and the arrangement of particles in space lead to the occurrence of a broad maximum in the scattering pattern (Figure 1, ROI V). We simulate these effects by implementing the adhesive hard-sphere structure factor, S_SHS(q) [35-37]. The interaction between particles at the distance x is approximated by the following potential U(x):

where R_HS is a hard sphere radius of particles, Δ is the width of a potential well, and τ is a stickiness parameter. The structure factor expression, S_SHS(q), is defined through the following set of equations (Equations 5–14). The structure factor is a function of four parameters in Equation 5, whereas the remaining variables and the associated expressions in Equations 6–14 are merely auxiliary by convention, so that we avoid writing a very long single expression. They should be consequently substituted into each other where suitable to obtain the final expression.

In this study, the stickiness τ is calculated from Equation 4, based on the literature data for silica–LZM composites synthesized under similar physicochemical conditions [24]: U(2R_HS < x < 2R_HS + Δ) = −2.5 k_BT; Δ = 0.1·(2R_HS) under the restriction that Δ_min ≥ 0.15 nm (i.e., the average H-bond length). The value of U in our study may slightly differ from literature values, because for constant pH value and salinity, the surface charge of NPs increases with decreasing size [38]. However, it is unlikely that U > −3k_BT, and within the considered range, the value of U will not affect our fitting results. ν is a local packing parameter, i.e., a local volume fraction within the aggregate, and for the random packing of polydisperse spheres it does not exceed 0.65 [39-42]. Figure 2 shows the scattering patterns at 0 and 100 s (Figure 1D) together with simulated curves based on Equation 2, in which the contributions of the silica NPs were taken from the Monte Carlo-fitted form factor (Figure 1B), while the contribution of S_eff(q,r_i) was included from Equation 5 (and Equations 6–14).

The simulations show the important effects that polydispersity has on the structure factor and the position of the correlation peak in ROI V. Typically, for correlations originating from (sticky) hard-sphere interactions, one considers the following dependence, for the approximated position of the peak at q_max:

Equation 15 infers that the expected average hard-sphere radius, R_HS, would be equal to the mean radius of a silica particle <r> ≈ 2.5 nm (Figure 1B, Table S2, Supporting Information File 1). This, in turn, suggests that silica NPs on average touch each other without any LZM molecules in between, or that the protein molecules, if present within the aggregates and among individual silica NPs, are very strongly deformed, likely to a point that they barely contribute to the determined R_HS. Nevertheless, the simulation in Figure 2A clearly shows that if the size distribution of silica NPs is actually correctly accounted for, then in order to fit the peak position accurately, an additional “spacer”, aR_eHS (additional effective hard sphere radius) has to be included in Equation 2 and Equations 5–14:

By setting merely r_i = R_HS,i (i.e., aR_eHS = 0) the position of the simulated peak visibly shifts towards higher q values with respect to the measured peak. Here, aR_eHS is associated with the presence of a single LZM molecule located in between individual silica NPs with the diameter of the LZM molecule represented by 2aR_eHS. The simulation in Figure 2B also shows that the packing factor ν within the aggregates, which directly correlates with the intensity of the broad peak around q ≈ 1 nm⁻¹, has to be relatively high (ν above ca. 0.4) in order to be able to simulate the intensity profile at q ≈ 1 nm⁻¹ in the later stages (ca. 100 s).

Low-q intensity increase, aggregation, and the structure factor expression of aggregates with internal correlations

The structure factor from Equations 5–14 does not reproduce the observed intensity increase at low q (Figure 1 and Figure 2A), because the sticky hard-sphere structure factor is derived under the assumption that the interactions extend to an infinite length scale, with respect to the probed volume. However, in our experiments, we initially have a finite number of “loose” silica NPs that are then rearranged to large aggregates upon mixing with lysozyme. Hence, microscopically such aggregates must have a finite size, even if their size (radii of gyration) cannot be determined directly from our scattering data, due to the used q range (see Experimental section). However, ultimately, to obtain a good fit an expression for S_eff(q) (Equation 2) has to account for both the fine structure within the aggregates causing the interparticle correlations (as in Equations 5–14) and also the low-q intensity increase due to the presence of the interface between the aggregates/composites and the solvent matrix.

A general expression for such an effective structure factor for particles within an aggregate/droplet was proposed originally by Hashimoto and co-workers [43]. Several variations and applications of this concept are furthermore known [37,44]. For the purpose of the analysis of our data, we further extended the expressions originating from Hashimoto et al. as we show below in a final form. In the Appendix section the complete derivation and the rationale are presented. Our derivation is essential here, because it allowed us to quantify indirectly the changes in the size of the aggregates, although the direct measurement of their radii of gyration was not possible. The so-derived general expression for S_eff(q,r_i) (Equation 17) is expressed as the sum between the structure function of an aggregate (“template”), S_agg(q), and the structure factor of the internal arrangement of the aggregate, S_int(q), which in our case becomes subsisted by S_SHS(q) (Equations 5–14):

where D is a fractal dimension describing the arrangement of primary particles within the composites, and A is a single collective fitting parameter in our model, which is proportional to the number density of aggregates, N_agg, and their specific surface area, SSA_agg. Hence, it expresses indirectly the size/extent of the aggregates. In Figure 3, we show that the fits with the introduced expressions for the partial structure-factor contributions indeed represent the structural features present in the selected scattering pattern (example @ 100 s).

However, as is evident from Figure 3A, the S_SHS(q) from Equations 6–14 has to be further improved, because in Figure 3A the correlation peak is relatively broad (“smeared out”), yet still intense. Typically one would expect such a broadened shape if the local volume-fraction parameter, ν, was smaller than derived from the best fit (i.e., below ca. 0.4). However, this would also inevitably yield a smaller relative intensity of this peak (see Figure 2B). Hence, in order to explain this contradiction, one has to remember that the position of the maximum and its shape are predominantly related to aR_eHS. The shape of the peak can be modelled substantially better if one allows for a distribution of this parameter in the fitting routine. The need for such mathematical treatment is in fact a manifestation of the actual physical effects, if we consider that aR_eHS represents a radius of a LZM molecule. A LZM molecule can become, at least partially, heterogeneously deformed (on average, in a global sense) within an aggregate, e.g., due to the variation of local forces, which, in turn, is a consequence of polydispersity of the silica NPs and the random character of the packing of the silica NPs. More importantly, since lysozyme is a small prolate ellipsoidal protein, with its principal semi-axes being 1.5 nm × 1.5 nm × 2.25 nm, the polydispersity in aR_eHS may account for the fact that the protein molecules can be differently orientated during adsorption to the NPs. Yet, so far we tried to represent their contribution through a (hard) spherical model. To overcome this, we used a Gaussian distribution to define the average structure factor <S_SHS,i(q)> in Equation 18. The application of this structure factor ultimately leads to smearing of the maximum at a constant value of ν and hence yields significantly improved fits (Figure 3B):

In Equation 18, the mean of the distribution (Equation 19) is <aR_eHS>, whereas σ denotes the standard deviation. This is the final expression used to represent S_int(q) in Equation 17 and to fit all the scattering curves from regions II and III in Figure 1. The numerical integration was performed for each i-th bin of the discrete size distribution characterizing the form factor (Figure 1B). The complete source code and selected data sets are deposited at GitHub.com [45].

Application of the scattering model

We used the above described scattering model to fit the time-resolved data set and to derive the fitting parameters as a function of the time (Figure 1, Period II, III and beginning of IV, 3010 curves). Figure 4 shows the time dependence of the five model parameters obtained by fitting the scattering data between 24.5 and 175 s.

During the first ca. 5 s after mixing, the values of <aR_eHS> and σ lie way outside reasonable error margins (Figure 4A,B). This is to be expected, because the corresponding local volume fraction, ν, is very low (well below 0.1, Figure 4C) at the beginning of the composite-formation process (i.e., the beginning of period II in Figure 1). Consequently, the contribution of <S_SHS(q)> to the structure factor during this period is mostly negligible with respect to S_agg(q). As ν reaches ca. 0.1 after about 30 s (transition between regions II and III in Figure 1) the actual evolution of <aR_eHS> and σ begin: <aR_eHS> starts at a value of ca. 0.7 nm and rapidly increases to ca. 0.9 nm after 50 s and then more gradually to ca. 1.1 nm after 150 s (period III in Figure 1). The associated standard deviation σ follows a similar trend as <aR_eHS>, growing from 0.35 to 0.5 nm after 50 s, and then levels off at ca. 0.5 nm after 150 s (end of period III in Figure 1A) within the fitting uncertainty. These initial rapid changes after up to 50 s are also reflected in the evolution of parameter ν (Figure 4C), where the parameter rapidly increases to ca. 0.45 and then remains constant within the fitting uncertainty after up to 150 s. Because <aR_eHS> directly relates to the size of lysozyme within the silica aggregate, we can link the changes in <aR_eHS> to possible changes in the structure/folding/shape of the protein. Between ca. 30 and ca. 50 s, where the LZM molecules rapidly induce the aggregation of the silica NPs, the LZM molecules appear to undergo a deformation (compression). However, as this process approaches equilibrium, through the internal densification of the aggregates, the molecules gradually return to their native dimensions. This result is in agreement with findings concerning the activity of lysozyme within composites with silica, showing that smaller silica NPs (as those used for our experiments) promote higher enzymatic activity of lysozyme, and that this, in turn, depends on the preservation of the native shape of the molecule upon composite formation [22].

Furthermore, the formation of the silica–LZM composite is dominated, from the very moment of mixing for the initial 20 s, by a rapid, 4.5-fold increase of parameter A (Figure 4D). This is best explained by the increasing number density of the aggregates, N_agg, and the associated increasing specific surface area, SSA_agg (see Equation 17 and Equation 29 in the Appendix). After t = 40 s, parameter A further increased, albeit at a slower rate. The concurrent evolution of the fractal dimension (parameter D; Figure 4E) suggests that initially (up to 50 s), the aggregates have a relatively open morphology with D < 2.4 and characterized by a limited contribution of <S_SHS(q)> due to ν < 0.1 (Figure 4C). Afterwards (after more than 50 s), the aggregates reached an internally denser state, as reflected by the steadying of both values for D (ca. 2.4) and ν (ca. 0.45). In other words, since these two parameters, D and v, reflect the internal structure of the aggregates from the perspective of the two structure-factor contributions (Equation 17), their evolution clearly indicates no further internal changes in the aggregates between 50 and 150 s. If such an internal densification processes had occurred, one would expect that it would have contributed to the decrease of the specific surface area of the aggregates, SSA_agg. Interestingly, however, parameter A (Figure 4D) keeps increasing after 50 s, i.e., after the internal dense structure is established, meaning that the product of number density of the aggregates and their specific surface area actually increased. This is possibly a result of secondary processes involving the “breakup” of larger aggregates into smaller units. Indeed, if we correlate the changes in A and ν (Figure 4F), we observe three stages of such secondary processes. In the first stage for ν < 0.1 (up to 25.5 s), A grows as a function of ν in a bound exponential mode, which translates into an increasing number of low-dimensional aggregates with hardly any internal correlations, forming an extended network of particles of low dimensionality D (Figure 4E). In the second stage, as ν increases from 0.1 to ca. 0.45 (25.5 to 70 s), A as a function of ν (Figure 4F) shows a linear dependence, indicating that as the number density of aggregates increases they also gradually densify, and that the growth of the aggregates occurs at the same rate as their internal densification. Finally, in the third stage, once ν remains relatively constant at ca. 0.45 (after 70 s), the product of the number density and the specific surface area of the aggregates continues to increase as documented by the increasing value of A, yet without any further dramatic changes to the internal structure/arrangement (i.e., constant values of D and ν), implying the aforementioned breakup of the larger aggregates into smaller units. These processes can be best explained as the initial rapid flocculation/clumping of NPs and LZM together into an extensive network just after mixing as the system is out of equilibrium, followed by the gradual evolution towards a steady state, in which smaller aggregates are more favourable.

During period IV (after more than 150 s), the time evolution of the three parameters (<aR_eHS>, ν and σ; Figure 4A–C) exhibited a characteristic discontinuity from the trends observed during periods II and III. This is because at times greater than 150 s, the scattering intensity at high q (which corresponds to the form factor, Figure 1D) changed significantly, so that the original form factor of pre-mixing silica NPs from Figure 1B was not representative for silica particles after 150 s. Thus, we could no longer use the fitted size distribution in our model, and any trends of these three parameters (Figure 4A–C) were not valid any more after 150 s. On the other hand, due to the fact that the low-q part of the data by definition is practically independent from the form factor, in fact the evolution of parameters A and D (trends in Figure 4D,E), even after 150 s are representative for the processes at the length scales corresponding to entire aggregates. However, due to the fact that in period IV our scattering model is no longer self-consistent, we did not analyse those trends.

Implications

The analysis of the evolution of the fitting parameters (Figure 4) draws a clear image of the four-step sequence of events during aggregation induced by the interaction between the protein LZM and amorphous silica NPs (Figure 5).

Immediately upon mixing, aggregation is induced due to the opposing surface charge of the silica NPs and the protein. An infinitely extensive and open (D = 1.8–2.2) aggregate network, from the point of view of the SAXS measurement, forms within ca. four data frames (ca. 200 ms). The so-formed network initially has no internal correlations, as is expected for a classical mass fractal [46]. However soon after (ca. 1 s), areas of correlated NPs–LZM domains start forming within the network and the increase in the internal volume fraction, ν, indicates an internal densification and ordering. This is also reflected by the fact that parameter D reaches a stable and relatively high value of ca. 2.4, which is characteristic for denser mass fractals. Such a fractal dimension for silica–lysozyme aggregates was previously reported [16,23] and can be associated with the diffusion-limited particle–cluster aggregation (DLPCA) mechanism [47-52]. The DLPCA growth mode is also evidenced not only by the plateau value of D, but also through the fact that the value of D increases with the size/extent of the aggregates [47,49] expressed indirectly by the parameter A (see region II in Figure 4D–E). This means that aggregates grow through the accretion of individual primary particles to larger aggregates [52], where aggregates as such become denser as their size increases, which in consequence favours the eventual occurrence of the correlated domains (which are the ultimate dense regions in the aggregate composed of smaller particles).

The parameters characterizing the interparticle correlation effects carry information about the size of the LZM molecules bridging the silica NPs. From Figure 4A–B it is clear that in this network the dimensions of the protein molecules are considerably smaller than the native dimensions of lysozyme in any possible orientation. Hence, this suggests that initially the binding of silica NPs by lysozyme involves a severe deformation/unfolding of the protein molecules, followed in time by a relaxation and increase in the protein dimensions towards a (more) native state (Figure 5). The final value of the radius of 1.1 ± 0.5 nm for the protein, which is reached before 150 s, is close to a radius of the protein in a side-on orientation (ca. 1.5 nm), rather than in the end-on orientation (ca. 2.25 nm). Our ex situ analysis of the dried composite samples (see Experimental section) showed that at 1000 ppm lysozyme, 32.7 wt % of the protein was incorporated into the composites. This means that for silica NPs precipitated at the concentration of 0.8737 g/L (calculated from the volume distribution in the SAXS patterns), the concentration of lysozyme in the composite was 0.4245 g/L. This is valid under the assumption that all available silica NPs were bound in aggregates with lysozyme. Hence, the number density of the protein molecules was N_LZM ≈ 1.8·10¹⁹ L⁻¹_, and for silica N_NP ≈ 1.8·10¹⁹ L⁻¹ (from SAXS by converting the volume distribution to a number distribution). This directly suggests that the silica–lysozyme aggregates are near almost stoichiometric, with one protein molecule associated with one silica NP. Such a stoichiometric relationship is actually expected for small silica particles of the size close to that of the protein molecule [22,53]. Su et al. [54] found that at small surface coverage the lysozyme attaches to silica NPs in a side on orientation, and recently the molecular dynamics simulations by Hildebrand et al. [55] also further confirmed that the side-on orientation of lysozyme with respect to silica constitutes the configuration of the highest attraction. This together with the relatively low dipole moment and the positive surface charge of the protein surface, potentially accounts for the bridging of the NPs by the LZM molecules, as the protein does not show a favoured orientation of the opposite active sites in the side-on orientation (i.e., both active sites show similar binding properties). In such a case one should indeed expect the DLPCA mode of aggregation, with the binding of the protein to the silica NP surfaces taking place through specific amino acids at the opposite sides of the molecule [24,55-57]. The densification of the internal structure of the aggregates reaches a steady point, when the LZM molecules relax to their native-like dimensions. Yet, at the same time the actual network constituting the composite, appears to break up into smaller aggregate units. The morphological changes of the composites further continue beyond 150 s. This is documented through the change in the form factor of the silica NPs that appear to grow in size, compared to the pre-mixed initial NPs. Although we cannot use our model to explain this last stage, we can speculate that the observed change is caused by a partial coalescence or Ostwald ripening of NPs inside of the aggregates. It is well documented that silica NPs synthesized from monosilicic acid are internally highly disordered and hydrated at their surfaces [31]. Previous studies found that the initially formed NPs if aggregated continued their growth, resulting in some cases in larger homogeneous silica spheres [31,58,59], and this may explain why silica is such a persistent scaling material in hydrothermal systems.

Conclusion

The analysis of in situ scattering data collected at 20 fps from the formation of silica–lysozyme composites showed that the processes included the formation of large aggregated structures in which individual silica NPs were bridged by LZM molecules. We developed and applied a new scattering model to underpin the changes of the morphology of the composites as a function of time. This model allowed us to unravel that the formation follows a diffusion limited particle–cluster aggregation (DLPCA) mechanism, which results in relatively densely packed mass-fractal-like aggregates within which non-fractal correlated domains of particles evolve. Furthermore, we used the scattering model to link the evolution in the measured structure factor to the lysozyme molecule, and we found that the aggregation processes involve severe deformation of the protein molecules, which is then followed by the relaxation towards the original dimensions.

Appendix

Derivation of Equation 17

The contribution of S_int(q) in Equation 17 has a rather simple form. However, the actual meaning of the parameter A extends beyond this short representation. Below, we present the derivation of the equation and we indicate the approximations we make on the way to the final form of the equation.

First we consider a structure factor function describing an arrangement of primary particles of a certain scattering length density (SLD), ρ_p, within an aggregate. The SLD of the space between the particles within the aggregate (the “template”) is ρ_t, while the SLD of the surrounding (“solvent”) is ρ_s. The average scattering contrast of the aggregate, Δρ_agg, is then:

where ν is a local volume fraction of particles within the aggregate (as in Equation 5). We also define any contrast fluctuation within an aggregate as:

If we write the structure function of an aggregate (“template”) as S_agg(q) and the internal arrangement within this aggregate as S_int(q), then similar to Hashimoto et al. [43] and Lin et al. [44], we can express the scattering intensity by the generalized Equation 22, in which “” denotes a convolution operation of the functions:

where V_agg and represent the volume and the volume fraction of aggregates in the solution, respectively. Here, S_int(q) = S_SHS(q) from Equations 5–14, whereas the S_agg(q) contributes to the increase in intensity at low q in the course of the formation of aggregates. The approximation in Equation 22 is valid when the overall radii of the aggregates, ξ, are significantly larger that the interparticle correlation distance from Equations 5–14 (ξ >> R_HS). Furthermore, we must re-normalize S_eff(q) in a different way than Hashimoto et al. and Lin et al. in order to combine it with the form factor defined in Equation 1 as for our scattering data the intensity for q→0 in a log–log representation does not level off to a finite value (Figure 2 and Figure 3). Consequently, it is impossible to determine V_agg. On the other hand the radii of the primary particles and, consequently, their volume, V_p, are known, and therefore it is reasonable to normalize S_eff(q) with respect to the primary particle volume rather than the aggregate volume:

Now let us assume that the aggregates follow mass-fractal behaviour and that the correlation function of mass fractals is described [60,61] by Equation 24:

where D is a fractal dimension, d is the Euclidean dimension (d = 3) and K is a normalization constant that is proportional to the mass and surface area of an aggregate. We do not include a cut-off function since in our case the intensity did not level off at low q as mentioned above. For the sake of simplicity let us assume again that the electron density of the aggregates is homogeneously distributed and its corresponding correlation function only involves a two-phase system. In this case the specific surface area of aggregates (SSA_agg) is proportional to the derivative of the correlation function at the near-zero length scale r→0:

This is valid for the considered small volume fractions (). The limit in Equation 25 exists essentially only for D = 4, and since D < 4, SSA_agg becomes increasingly larger for decreasing length scales, yielding an infinite surface area at an infinitely small length scale r. However, since the aggregates are composed of primary particles with a typical radius, R_HS, we can say that the aggregate does not contain smaller features than those primary particles (i.e., r ≥ R_HS). Therefore, we find a finite specific surface area for mass fractal aggregates for lim(r→R_HS):

This way this newly derived K constant in Equation 26 substitutes the K constant from Equation 24 and thus, we can use the correlation function from Equation 24 to calculate the structure factor. Please note the Hashimoto et al. described the structure factor in such a way that it is normalized as a form factor, i.e., it is normalized by the total volume. This is in line with the structure factor of mass-fractal aggregates as described by Sorensen and Wang [60], yet it is different from a better-known derivation by Teixeira [62]. Both Sorensen and Wang’s and Teixeira’s approaches are valid as long as one considers normalizations explicitly. The structure factor is described by the rotation-averaged Fourier transform:

where N_agg is the number density of aggregates and SSA_agg is their specific surface area. By substitution of Equation 27 into Equation 23, we obtain:

For spherical primary particles V_p is known and a final form of the equation for the effective structure factor could be derived:

In Equation 29 we introduced several simplifications. Firstly, N_agg, SSA_agg, Δη and Δρ_agg are essentially unknown, and it is impossible to determine any one of these parameters independently. They have to be combined into a collective parameter. This is necessary as Equation 29 was derived for a system characterized by a monodisperse particle distribution with only a single value of for a given R_HS and the resulting N_agg. For a polydisperse distribution, as in our study, although a population of primary particles is described by D(r_i) (Equation 2), the resulting distribution of aggregate sizes will be totally independent from this initial distribution, and it will also be unpredictable. Secondly, the R_HS^1−D component in Equation 29 could be important, since it determines the high-q cut-off at which the contribution of the structure to the intensity lessens, and the form factor dominates. In Equation 29 this very transition point is dominated by a S_int(q) contribution and its strong correlation peak. Hence, we assume that R_HS^1−D ≈ 1. Thirdly, the remaining part of the expression depending on parameter D is practically constant at a value of ca. 1.2, and although we could introduce it explicitly in the model it does not affect the final trends. Hence, as a result of the above approximations, we use A as a single collective fitting parameter in our model. Changes in A therefore should be interpreted primarily as the average change of the product of the number density of aggregates and their specific surface area, and these two physical parameters are related to the size (or “extent”) of the aggregates.

Experimental

Synthesis of amorphous silica–lysozyme composites

Sodium metasilicate (Na₂SiO₃·5H₂O, technical grade), hen egg-white lysozyme (crystalline, powdered, >90% pure, residual sodium acetate and chloride) and HCl (37%, analytical grade) were purchased from Sigma-Aldrich. Separate stock solutions of dissolved silica (SiO₂ = 1000 ppm, pH 12.5) and lysozyme (5 wt %, pH 3.5) were prepared by dissolving the required amount of sodium metasilicate or lysozyme in ultrapure deionised water (18.2 MΩ·cm). Silica NPs were prepared in a 500 mL plastic reactor by neutralizing the silica stock solution through titration with HCl until pH 7.5 was reached. This neutralized solution was left to polymerize and age for 16 h. Silica–protein composites were obtained by mixing the silica NP solution with a pre-measured amount of the lysozyme stock solution under rapid stirring (500 rpm) to yield a SiO₂ NPs solution with 1000 ppm lysozyme (final pH 6.9, salinity 20 mM).

Scattering experiments

The formation process and the development of the structure of the silica–LZM composites was studied in situ and in a time-resolved manner by using synchrotron-based small-angle X-ray scattering (SAXS) at the BioSAXS beamline [63] P12 of the EMBL at PETRA III (DESY, Germany) using a monochromatic X-ray beam at 10 keV. Two-dimensional scattered intensities were collected at small angles with a Dectris Pilatus 2M (2D large-area pixel array detector) using an acquisition time of 50 ms per frame. Transmission was measured by means of a photodiode installed in the beam stop of the SAXS detector. A sample-to-detector distance of ca. 3 m allowed for a usable q-range of 0.04 < q < 4.5 nm⁻¹. The scattering-range at small angles was calibrated against silver behenate, and the intensity was calibrated to absolute units against water. For the in situ experiment, first, the starting silica NP solution was continuously circulated between the reactor (where the suspension was stirred at 500 rpm) and the flow-through cell with embedded quartz capillary (ID 1.7 mm, wall thickness 50 µm; aligned with the X-ray beam) using a peristaltic pump (Gilson MiniPuls 3, flow ca. 500 mL/min; tubing: ID 2 mm, total length 2 m; reactor-to-cell-distance: 0.7 m of tubing). All experiments were conducted at 21 °C. Once a SAXS baseline for the silica NP solution was recorded, the pre-measured amount of lysozyme stock solution was pumped into the reactor at a fast rate. This injection was done remotely from the operator hutch via a 10 m long PTFE tube (ID 4 mm) that was routed into the reactor located in the experimental hutch. The tube was filled in such a way that the LZM solution was located in the last ca. 40 cm of the tube on the reactor side. The other end of the tube in the operator room was equipped with a 50 mL syringe filled with air. Thus the experiment started with recording of 24 s (480 × 50 ms) SAXS patterns of the silica NP solution circulating through the capillary, prior to the fast injection of the entire content of the tube containing the lysozyme with a single rapid push of the syringe plunger that lasted ca. 200–400 ms. This fast injection rate in combination with the fast stirring in the reactor (500 rpm), pumping (500 mL/min) and fast data acquisition (50 ms per SAXS pattern) provided the best possible conditions for the characterisation of all the steps leading to the formation of the silica–LZM composites. The used experimental set up introduced an unavoidable dead-time of ca. 500 ms between the injection moment and the first actual measurement of the mixed solution, i.e., the time required for the mixed solution to reach the capillary where the SAXS pattern was recorded. In order to be able to analyse and model the silica–LZM composite scattering patterns we also acquired a series of backgrounds and reference samples including an empty capillary and a capillary filled with water, silica stock solution, LZM solutions at different concentrations. The initial SAXS data processing and reduction included a series of automatic post-data-collection steps including masking of undesired pixels, normalizations and correction for transmission, instrumental background subtraction and data integration of the collected 2D data to 1D. Further data processing and water background subtraction, model fitting, validation and analysis, were performed through a custom-made script developed in GNU Octave [64,65]. The script we developed as well as all the documentation and the selected scattering curves are available at [45]: https://github.com/tomaszstawski/SilicaLysozymeSAXS. In a first instance for the model, we obtained the size distribution of the initial silica NPs from a Monte Carlo fitting implemented [29,30] in MCSAS under the assumption that the silica NPs particles were spherical [31].

Characterization of dry samples

To cross-correlate the in situ SAXS data, the silica–lysozyme suspensions were dried in an oven at 40 °C for ca. 48 h. The resulting powders were washed five times with MilliQ water to remove excess lysozyme and salts followed by a 2nd drying step at 40 °C. The amount of lysozyme associated with the composites was quantified by determining the total carbon content in solids by mass spectrometry (DELTAplusXL ThermoFisher) with a Carlo-Erba NC2500. From these analyses the lysozyme content was calculated using the molecular formula C₆₁₃H₉₅₉N₁₉₃O₁₈₅S₁₀ and the molecular weight of 14313 g/mol for lysozyme [66] (ProtParam based on UniProtKB entry P00698).