{'fraction_non_alphanumeric': 0.023133462282398452, 'fraction_numerical': 0.005609284332688588, 'mean_word_length': 4.69143549097314, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A hard-sphere model of protein corona formation on spherical and cylindrical nanoparticles

A nanoparticle (NP) immersed in biological media rapidly forms a corona of adsorbed proteins, which later controls the eventual fate of the particle and the route through which adverse outcomes may occur. The composition and timescale for the formation of this corona are both highly dependent on both the NP and its environment. The deposition of proteins on the surface of the NP is related to processes of random sequential adsorption and, based on this model, we develop a rate-equation treatment for the formation of a corona represented by hard spheres on spherical and cylindrical NPs. We find that the geometry of the NP significantly alters the composition of the corona through a process independent of the rate constants assumed for adsorption and desorption of proteins, with the radius and shape of the NP both influencing the corona. Moreover, we demonstrate that in the condition of strong binding such that the adsorption is effectively irreversible the corona content reflects the protein mobility and concentration in solution rather than their binding affinity.


INTRODUCTION
Matter, even elements which are well-known on everyday size scales, takes on new properties on the nanoscale due to increased surface area and quantum mechanical effects, such that even a famously inert metal like gold may become highly reactive in the form of nanoparticles (NP) (1). These new properties, however, potentially lead also to hazards when the NPs are inhaled, absorbed through the skin, ingested, or otherwise find their way into the human body (2,3). A key finding so far has been the observation that a NP in a biological medium rapidly forms a "corona" of adsorbed proteins and other molecules (4,5). This corona is known to be linked to the transport of the NPs around the body, their uptake by cells, and potential adverse outcomes (6)(7)(8). Thus, the prediction of the content of this corona is key to predicting the safety of novel nanomaterials.
Qualitatively, the corona may be separated into two classes of bound protein: "hard" and "soft" (4). The hard corona refers to the proteins strongly bound directly to the surface of the NP, whereas the soft corona is a collection of proteins which are loosely bound on top of the adsorbed protein layer. The content of the hard corona is highly dependent on the material and morphology of the NP, and evolves through the adsorption of proteins upon their diffusion from the bulk, followed by their later desorption and replacement by other proteins over several hours (4,5). The rate at which proteins desorb is linked to the strength with which they are bound to the NP. The weakly bound proteins are typically smaller and so diffuse faster through the biological medium. Thus, they have a higher rate of collisions with the NP and may feature more in the earlier stages of the corona before being replaced by more strongly-binding proteins which collide less frequently, but do not desorb as rapidly. The phenomenon of protein replacement in the adsorption layer is known as the Vroman effect (9). The chemical composition of the NP and any modifications to the surface directly influence the strength with which different proteins bind to it, for example, titania NPs are known to be strongly hydrophilic and thus bind more strongly to certain classes of proteins (10). Intuitively, the curvature of the surface of the NP also influences the strength of the binding: a flat NP will have more of its surface contacting a Manuscript submitted to Biophysical Journal 1 arXiv:2104.13699v1 [cond-mat.soft] 28 Apr 2021 given protein, and thus experience a stronger attractive force than a curved NP that makes a point contact. We may readily expect therefore that the strength of the binding differs between NPs of the same material but different surface curvatures, i.e., as a function of the geometry of the NP, as is observed in calculations of the binding energy through coarse-grained methods (11,12). Beyond this, we expect that the curvature of the surface plays a role in defining how tightly proteins may pack on the surface as a purely geometrical effect (13,14). Indeed, as observed in recent work (15), spherical and cylindrical NPs of the same material may have have differing corona compositions, which are attributed to a combination of these two effects.
Given the relevance of the corona to biological outcomes, it is therefore essential to get insight into the mechanisms of its formation and predict its content. This can be done experimentally through a variety of methods, but these are typically highly sensitive to the exact concentration of proteins used, which may vary significantly between the samples of real biological fluids, and cannot be immediately applied to novel NPs. We therefore turn to theoretical and computational methods to predict both the evolution and the steady state. In principle, this could be achieved by performing a molecular dynamics simulation of the NP immersed in the biological liquid. A brute force solution of the problem is difficult due to the combination of extremely large numbers of atoms and long timescales involved. Atomistic simulations of individual proteins typically require several days to produce a few nanoseconds of simulated time, while the corona may consist of hundreds of proteins and develop over several hours to days (16). Thus, to model the actual adsorption-desorption kinetics it is necessary to employ a coarse-grained approach in which the biomolecules are represented by less complex structures that are more amenable to simulation. Even these coarse-grained simulations may be still limited in terms of the amount of computational time taken, and so to probe the dynamics of corona formation on biological timescales one may need to turn to alternate methods, e.g. rate-equation treatments or different simulation techniques (17)(18)(19)(20)(21)(22)(23)(24)(25).
Previous work on the formation of the corona based on rate equation treatments has primarily focused on a mean-field (MF) approach in which the rate of incoming proteins is reduced by a factor dependent on the total surface area occupied by proteins but further steric factors are not taken into account. That is, if 50% of the surface area of the NP is occupied by a protein, the binding rate for all incoming proteins is reduced by 50%. This assumption is justified if the proteins are sufficiently deformable that they pack perfectly tightly without any gaps. If this is not the case, then it is straightforward to envision a situation in which the total uncovered area is sufficiently large to allow the binding of further proteins, but is distributed across a number of gaps each of which are too small individually to admit a further protein. In general, proteins cannot arbitrarily adjust their bound configurations and so it is reasonable to assume that some parts of the NP will be uncovered but unavailable for binding. If we consider in particular globular proteins which may be approximated as spherical, then it is clear that there must be gaps between these proteins due to the impossibility to perfectly pack spheres. This leads us to consider a model of random sequential adsorption (RSA) with the addition of desorption and surface diffusion, see e.g. Ref (26) for an overview of models. In this paradigm, the proteins are modelled as hard spheres with well-defined positions on the NP. Intuitively, we may expect this hard-sphere (HS) model to produce different results to the MF approach. Indeed, numerical simulations have shown that the HS model predicts effects arising purely from the geometry of the NP and deviations from the mean-field predictions (13). Provided that the rate equation method can be adjusted to include these effects, in principle this enables a more physically realistic modelling of the corona formation. However, rate equation treatments have a limitation in that the medium in which the NP is immersed may contain a very large number of different protein species, each of which may bind at an essentially arbitrary orientation relative to the NP with an affinity unique to each orientation. Numerical integration of such a large number of dependent variables rapidly becomes infeasible, while models more complex than the mean-field typically lack analytical solutions for even a single spherically symmetric protein. Thus, we must consider also computational approaches capable of simulating large periods of time, e.g., kinetic Monte Carlo approaches in which the time between the adsorption and desorption of proteins is not directly simulated.
In this work, we employ a combined analytical and computational approach to investigate the effects of the geometry of the NP on the corona composition in terms of how this alters the packing efficiency and the rate at which proteins adsorb to the surface. We demonstrate that the mean-field and hard sphere (HS) models predict significant differences in the timescales and numbers of adsorbed proteins for a simple system of three common blood serum proteins previously employed in studies of corona formation on spherical NPs (17,18). Using an approach based on scaled particle theory as previously proposed for adsorption onto planar surfaces (27), we develop a rate-equation treatment for the HS model that may be numerically integrated to obtain the evolution of the corona on spherical and cylindrical NPs, finding that the geometry of the NP alters the composition of the corona even when the change in binding affinity due to the surface curvature is neglected. We confirm these results using Kinetic Monte Carlo simulations, finding excellent agreement. Our model highlights the importance of geometric effects, binding strength, and the surface diffusion of the bound proteins.

Rate equations
We consider a model of reversible adsorption of proteins to the surface of an NP in which the protein is in physical contact with the surface and occupies an area thus making it inaccessible for other proteins. The surface area of the NP is NP , and we denote the number of available binding sites for a given protein as = NP / . We define the number of proteins of type which are bound to the NP to be and further define the surface coverage of this type to be = / . Each protein is further characterised by a rate constant for adsorption, , (in units s −1 M −1 ) and desorption , and a bulk molar concentration [ ].
We assume that this bulk concentration remains fixed during the evolution of the corona due to either diffusion of proteins from a reservoir or a sufficiently low concentration of NPs that the effect of binding on the concentration of unbound proteins is negligible.
Once a NP is immersed in a protein solution, the evolution of the resulting corona using a rate-equation treatment is given by where proteins of type collide with the NP at a rate given by , [ ] and (N) is the probability that a colliding protein of this type successfully binds to the surface of the NP as a function of the set of bound proteins, denoted by N = ( 1 , 2 ... ).
We further assume that the protein binds to the surface with unit probability if there is sufficient space for it to do so. Thus, the probability of binding is determined by the probability that a protein randomly positioned on the surface of the NP does not overlap with other proteins. In the mean-field approach of Refs. (17,18), this probability is simply equal to the fraction of the surface area of the NP that is not already covered by proteins, This function is independent of the identity of protein and thus the rate of binding of all proteins is reduced by the same amount. Using this definition and changing variables from to we obtain which is independent of the set of values . Thus, in the mean-field model, the corona composition defined in terms of the surface coverages depends only on the concentration of unbound proteins and the rate constants for adsorption and desorption. The set of linear first-order differential equations given by (3) can be solved analytically using standard techniques. Here, we simply give the steady-state corona composition : If the proteins are reasonably strongly binding or have a high concentration then the occupied fraction of the NP surface rapidly approaches unity, that is, the NP becomes completely covered in proteins with no gaps between them. This simple model is valid if the area occupied by a protein can be deformed to enable this optimal packing. In general, however, proteins possess some degree of rigidity preventing this deformation. Consequently, gaps will likely exist between proteins bound to the surface.
In the mean field model, the sum of these gaps represents valid area for the binding of further proteins, even if each gap is individually too small to admit a protein. A further limitation of this model is that it assumes that the coverage of the NP reduces the binding rate of each protein identically. Intuitively, we might expect that small proteins are less effected by the increased coverage due to their ability to fit better into gaps between proteins already present on the surface. To quantify these effects, we investigate a model in which the projections of the proteins are taken to be rigid and the proteins are modelled as spherical. We refer to this as the hard-sphere (HS) model, although strictly speaking it is the projections of the proteins which are hard in the sense they cannot overlap or deform. As we are interested only in the monolayer hard corona, the binding process is equivalent to the two-dimensional random sequential adsorption and desorption model which has been thoroughly investigated for planar surfaces. In particular, it has been shown that the acceptance probability for an incoming particle in this model can be approximated by (27) Unlike in Ref. (18), this is the exact solution to the steady-state as we have made the approximation that the concentration of unbound proteins remains fixed.

Manuscript submitted to Biophysical Journal 3
where is the excess chemical potential for the insertion of a new protein of type as a function of the current state, under the assumption that the state is currently in equilibrium. As before, a protein successfully binds to the NP with probability and rejected with probability 1 − . In the case of the adsorption of disks onto a planar surface, it was shown that this assumption is generally quite good even for non-equilibrium states and that the required acceptance probability can be obtained from scaled particle theory (SPT) (27). The probability for a freely-rotating convex particle to have sufficient space to bind is given by where and are the area and perimeter of the projection of protein onto the surface, s = ( 1 , 2 ... ). This expression implicitly assumes that the bound particles diffuse on the surface, and we investigate the effects of this assumption later.
In general, the acceptance probability in the HS model is reduced in comparison to that in the MF model, reflecting the fact that an incoming protein is not necessarily able to find a suitably large gap to bind to the NP even if there is sufficient space in total. To illustrate this effect, we consider the special case of a single type of spherically symmetric protein binding to a spherical NP, for which we can assume the projection of the molecule is circular. In this case, Eq.(6) becomes independent of the actual sizes of the protein and NP, but does not reduce to the MF case. In Figure 1 we plot the average time taken to adsorb a new protein, 1/( [ ] ) as a function of the surface coverage for both the mean-field and hard-sphere cases, setting [ ] = 1 −1 such that one potential adsorption occurs per second. It can immediately be seen that the time taken to adsorb a single new protein rapidly increases even at modest surface coverages in the HS model compared to the mean-field model. Since the formation of the corona likely involves many adsorption events, it follows that even if individual adsorption attempts proceed quicker than the rate of 1s −1 considered in this simple example, the overall formation of the corona will proceed slowly once a moderate surface coverage is reached. Consequently, the steady-state may not be reached on experimental timescales unless desorption of bound proteins is very fast and the equilibrium surface coverage is low. A rough estimate of the total collision rate coll = [ ] for a spherical protein to a spherical NP can be obtained from diffusion theory (28), where is the diffusion coefficient. For a typical protein diffusion coefficient on the order 1 m 2 s −1 (and neglecting the diffusion of the larger NP), protein concentration 1 M, with + NP = 20 nm, this implies a few hundred proteins collide with a given NP per second. With this collision rate, we may reasonably expect surface coverages on the order of 0.7 in the timeframe of a few hours.  It is reasonable to assume that we can employ Eq.(6) to estimate the probability for the successful insertion of proteins into the corona of curved NPs, provided that we correctly calculate , to reflect the projected area of the protein onto the surface of the NP when evaluating the insertion probability. Under this assumption, the evolution of the corona in the hard-sphere model is given by where will in general be a function of both the geometry of the NP and of the proteins. We may therefore expect that, unlike the MF case, the corona composition in this model will depend on the geometry of the NP even if the parameters , , [ ] are kept fixed, since this will alter the values of , used to calculate the value of . For a single protein, we can calculate the steady-state surface coverage as a function of the equilibrium constant and concentration, as shown in Fig. 1 and compared to the equivalent prediction for the MF model. It is immediately apparent that this steady-state coverage is much less in the HS model than in the MF, with two key implications. Firstly, we can generally expect fewer proteins to be bound in the HS model, even at very high concentrations or if the binding is particularly strong. Secondly, given an observed surface coverage (e.g. obtained from experiment) and initial concentration of proteins, different equilibrium constants would be obtained depending on which of the two models was employed, as noted in Ref (13).

Binding areas
To apply this model, we must define the area occupied by a protein on the surface of an NP. We consider two morphologies of NP: spherical NPs of radius NP and rod-like cylindrical NPs of radius NP and effectively infinite length. Planar surfaces can be obtained as a limit of either of these with NP → ∞. We assume that all proteins are approximately spherical and characterised by a radius , and that a bound protein is in physical contact with the NP. We project the protein onto a curved surface as shown schematically in Fig. 2. For a spherical NP, this occupied area is given in terms of the angle = sin −1 such that For a very large NP, the area tends towards 2 , i.e., the maximum cross-section of the protein. Conversely, for an NP much smaller than the protein, the number of available binding sites approaches the limiting value of 2, in agreement with the geometrical constraint that only two spheres may touch at a single point. In all cases, the number of binding sites is given by the surface area of the sphere divided by the area occupied by a protein, = 4 2 NP / . For the HS model we also require the projected perimeter of the protein. By reasons of symmetry, it follows that the projection should be circular and thus the perimeter can be obtained from the area by = 2 √︁ / . Next we turn to a cylindrical NP. The irreversible adsorption of monodisperse spheres onto a cylinder has been recently investigated in Ref. (14) and demonstrated to be equivalent to the deposition of approximately ellipsoidal shapes onto a planar surface, as depicted in Fig. 3. The exact area of the projection of the protein onto the cylinder is given by integrating the expression in Ref. (14) and re-scaling to restore the dependence on NP , where ( ) and ( ) are the complete elliptic functions of the second and first kind respectively and are given in terms of the parameter = 2 . The area can be used to calculate the number of available binding sites for a cylinder of finite (but large) length and convert the surface fraction into the number of bound proteins according to = 2 NP / . To apply the HS model, we must take into account the fact that the projections are not freely-rotating and thus Eq. (6) does not apply. Following Refs. (29,30), the acceptance probability for shapes with a fixed orientation is given by Manuscript submitted to Biophysical Journal 5 Figure 2: The projection of a spherical protein (blue) of radius onto an NP (gray) of radius . The angle is defined as the angle between the line connecting the centres of the NP and protein and the line tangent to the protein and passing through the centre of the NP. The area shaded in pale blue is assumed to be occupied by this protein and unavailable for binding by other proteins, blocking off a region of the surface of the NP. Under this projection, spherical proteins adsorbed to the NP as shown in the centre image can be mapped to two-dimensional regions on the surface of the NP as shown on the right.
is the excess excluded area for protein due to protein , defined in terms of the co-area which is the area of space in which the centre of protein cannot be placed without causing an overlap with protein . For the particular case of binding to a cylinder, we approximate this area by the area of an ellipse with axes + ,˜+˜, wherẽ and = ( + ) (˜+˜). It is important to note that, although we have presented the binding areas assuming a spherical protein, the general model proposed here depends only on the projection of the protein onto the surface and not its exact shape. Since Eq.(6) is applicable to freely-rotating convex projections and Eq.(12) to projections of fixed orientation, a wide range of proteins can be modelled using this approach provided that their projected area and perimeter can be calculated. Thus, a complex protein can be treated by estimating binding constants and an associated binding area for each possible orientation relative to the surface of the NP. In this work, we limit ourselves solely to considering spherical proteins to simplify the computational approach described later. We also note that it is the projections of the proteins which cannot overlap, which is a more strict requirement than simply that the proteins themselves do not overlap. This assumption is necessary both to allow small proteins to block larger ones (31) and to ensure that small proteins cannot fit in the gap between larger NPs and the surface of the NP. Physically, this assumption represents that the protein makes contact with the NP over a surface greater than a single point.

Computational methods
The set of equations in the MF approximation can be integrated analytically to produce the evolution of the corona as a function of concentrations and shapes of the proteins. The more complex HS, however, must be integrated numerically. We employ the numerical integration routines in Mathematica 12.1 to do so (32). To validate these equations, we also perform simulations of the adsorption process. As discussed previously, even highly coarse-grained models of corona formation using molecular dynamics take an extreme amount of computational time. Thus, we instead employ a kinetic Monte Carlo (KMC) method (33) that can readily be extended to large numbers of proteins. This algorithm samples a set of events -adsorption or desorption of proteins -advancing from one to the next without requiring the evaluation of the time inbetween them, allowing for the efficient simulation of a large number of events and thus an extended amount of time. At a given time , the algorithm generates a list of possible events, randomly selects one with a probability weighted by the rate at which that event occurs and the system is updated accordingly. Afterwards, the time is advanced by a random amount given by where is the sum of all event rates and is a random number uniformly distributed in the interval [0, 1]. This process iterates until the time reaches a pre-determined stopping point. The possible events correspond to proteins colliding with (and potentially adsorbing to) or desorbing from the surface of the NP. Collisions occur with a rate given by , [ ] and desorption with a rate , . If a desorption event rate is selected, a randomly chosen protein of the specified type is removed from the surface of the NP. Conversely, if a collision occurs, the protein is either accepted or rejected with criteria depending on the model employed. In the MF model, the current occupied surface fraction is calculated, and a random number in the interval [0, 1] is generated. If this random number is greater than , the protein is accepted and added to the state of the system. If not, the simulation proceeds to the next event. In the HS model, the location of each protein is explicitly tracked by generating a location for each colliding protein, where these positions are uniformly distributed on the surface of the NP. In the spherical case, a pair of angles , are randomly generated to produce a uniform distribution of points on the surface of the sphere. This is achieved by drawing two random numbers 1 , 2 from the uniform distribution [0, 1] and taking = 2 1 , = cos −1 (2 2 − 1) (34). A protein of the corresponding type is then inserted at the Cartesian coordinates given by ( NP + ) [cos sin , sin sin , cos ] provided that it does not overlap with any existing proteins. An overlap is deemed to have occurred if the condition, cos −1 (sin sin + cos cos cos( − )) < sin −1 NP + is found to hold for any protein already adsorbed to the NP and where protein is the new protein. This ensures that the projections of the proteins onto the surface do not overlap, as assumed in the rate-equation model, and is derived by projection of the centre of the two proteins onto the sphere and calculating the minimum angular distance between them allowable without overlap occuring.
For the cylindrical NP, the position of the protein is defined by the coordinates , , , where = NP + such that the protein is in contact with the NP and binding occurs only to the curved surface. The angular coordinate is chosen from a uniform distribution in [0, 2 ]. The final coordinate is chosen uniformly from the interval [− /2, /2], where is the length of the cylinder. To minimise edge effects, we ensure that is large compared to the typical radius of the proteins and employ periodic boundary conditions at ± /2. Unlike in the spherical case, we cannot write down a simple expression to test if a pair of proteins overlap due to the more complex shape of the projection of proteins onto the cylinder. Instead, for a pair of proteins of radii , , we determine if overlap occurs by setting the radial co-ordinate for both to , = NP + max( , ) and determining if this causes a physical overlap between the proteins, i.e. if the distance between their centres is less than + . This procedure ensures that a smaller protein cannot bind in the region between a larger protein and the NP, which would be possible if physical overlap was only checked for at their actual location, is consistent with the definition used for the spherical NP, and avoids a computationally expensive test for overlapping based on the projections of the two proteins onto the surface. By projecting the locations of proteins obtained from simulations onto a planar surface as depicted in Fig. 3 we have determined that this test correctly prevents the overlap of the projected areas of proteins for those considered here, although in general this method may fail if the proteins are much larger than the NP. To implement the periodic boundary conditions with respect to the ends of the cylinders, overlap is tested for the incoming protein against the original set of bound proteins and copies of this set translated by ± along the axis of the cylinder. For both spherical and cylindrical NPs, these overlap criteria are expected to produce equivalent results to the theoretical expressions Eq. (6) and Eq. (12) respectively, but the KMC program does not evaluate these pobabilities when assessing if a protein is accepted or not. Thus, this method provides a means to test the validity of these expressions for the systems considered here.
For the majority of the simulations, we assume that once a protein has adsorbed to the NP it remains fixed in place until it desorbs. It is known, however, that the diffusion of proteins on the surface of the NP increases the packing efficiency (35), and the rate-equation model assumes that the adsorbed proteins are effectively in a fluid phase (27). The desorption and readsorption of proteins produces a restructuring of the surface equivalent to the diffusion, justifying the exclusion of this process. If, however, proteins are sufficiently tightly bound that they do not escape the NP but are still free to diffuse across the surface, then it is of interest to investigate how the surface restructuring impacts the results. To implement this, after each event the proteins are selected in a random order and the angular position of each molecule perturbed by a small random amount. These perturbations are randomly drawn from zero-mean normal distributions with standard deviations given by (36) where Δ is the time between the most recently applied event and the event before and is the surface diffusion coefficient, which we set to 10 5 nm 2 ·s −1 for all proteins. A test move is generated and accepted if it does not result in an overlap between proteins and otherwise rejected, with up to four test moves attempted before moving to the next protein, with proteins selected in a random order. This requires re-checking for collisions following each trial perturbation and so dramatically slows the simulation. Thus, we do not employ this procedure for most simulations and have implemented it in the code only for spherical NPs.

RESULTS
In this work, we are primarily interested in determining how the geometry of the NP influences the corona, independently of the strength with which a given protein binds to the surface of the NP. Thus, we vary the size and shape of the NP while keeping the concentrations and rate constants for adsorption and desorption for each protein fixed and present results comparing the KMC simulations to the rate-equation model. For consistency with the literature, we primarily employ the same system of proteins as in Ref. (17) using the same rate constants and concentrations quoted there as summarised in Table 1. This system consists of three proteins: HDL, HSA, and fibrinogen (Fib). Of these, HDL and HSA are globular and so can be reasonably modelled as spherical. Fibrinogen is a rod-like protein, which can be treated using the analytical model but not in the current implementation of the CoronaKMC program for side-on orientations. For the purposes of this work, we simply assume the binding is achieved through an end-on configuration, such that the binding profile is approximately equivalent to that of a globular protein. The exact identity of these proteins is immaterial as we are focused on the mathematical treatment rather than attempting to match experimental results, and we simply require a range of concentrations and binding constants. Thus, Fib can also be considered to be a large but weakly binding globular protein, with the name used purely to ensure consistency with Ref. (17). Of the proteins considered here, HDL binds significantly more strongly than the other two, with a characteristic lifetime 1/ on the  Table 1: Summary of the parameters for the three proteins HDL, HSA and Fib used for the numerical investigation taken from Ref (17). The proteins FP1 and FP2 are fictional proteins based on modifications of the original three proteins to highlight the effects of varying parameters in the model. order of ≈ 9 hours compared to ≈ 8 minutes for HSA and Fib. In certain simulations, we also consider two fictional proteins corresponding to one comparable to Fib but with a higher concentration, and one which is a slightly larger form of HDL. To enable a comparison between NPs of different sizes and morphologies, we present all results in terms of the surface coverages , using the relation = / to convert the numbers of bound proteins observed in KMC simulations to the surface coverage. First, we confirm that the KMC algorithm can be used to generate results in the mean-field approximation by performing simulations using this algorithm and comparing it to the analytical solution to the rate-equation treatment. We find that the KMC method produces results essentially identical to the analytical method, as can be seen from the clustering of the trajectories around the analytical prediction for spherical NPs of radius both 5 and 50 nm in Fig. 4. Moreover, it can be seen that the mean surface fraction is independent of the radius of the NP, although the fluctuations around this number are larger for the smaller NP. The coverage for single-particle adsorption on infinite surfaces is known to follow a Poisson distribution (37), for which √︃ 2 − 2 ∝ √ , such that fluctuations in the surface coverage scale as √︁ / . The true distribution for multi-component adsorption on finite surfaces is of course likely to be much more complex, but heuristically we can expect similar scaling of fluctuations, especially when the steady-state corona consists of essentially only a single protein as in the present case.  Fig. 5. It is immediately apparent that this significantly alters both the surface coverage and the rate at which the system approaches equilibrium compared to the mean-field model. The analytical model using the surface coverage given by Eq. 6 is in good agreement with the numerical results. It can be seen that, unlike in the MF model, the surface coverage in the steady-state depends on the radius of the NP. Since we have left the rate constants unchanged, this implies that this is a purely geometric effect arising solely from the different packing efficiencies of spheres on curved surfaces. Likewise, in Fig. 6 we present a comparison between the numerical simulations and the numerically-integrated rate equation treatment for cylinders of radius = 1 and = 5nm, again finding that the radius of the NP plays a role in defining the evolution of the corona even when rate constants are not altered. In both cases, the time taken to reach the steady-state is on the order of several hours despite the relatively quick desorption of individual HSA molecules, which occurs on the order of minutes. This is most likely caused by the preferential re-adsorption of HSA into the gaps left by its desorption compared to the adsorption of the larger HDL molecules. Certain nanomaterials, for example gold, are known to adsorb proteins so strongly that the binding is essentially irreversible, , ≈ 0 (12). In this case, we may intuitively expect that the corona primarily represents the proteins which collide with the NP more frequently than others. The rate-equation model proposed previously is in principle still valid if we set = 0 for all proteins, however, the probability for successful adsorption of a protein in the hard-sphere model assumes a fluid-like phase of bound proteins (27). If proteins cannot desorb, then the bound layer of proteins does not meet this requirement and so we may expect deviations between the simulations and rate-equation model. The fluid-like behaviour can also be obtained by allowing bound proteins to diffuse on the surface of the NP and so we perform simulations including the effects of surface diffusion. The results in Figure 7 show a comparison between the simulations with (points) and without (open circles) surface diffusion for the protein set in the absence of desorption. Clearly, the rate-equation model remains accurate when surface diffusion is enabled for irreversible adsorption, but fails to accurately describe the results if the adsorbed proteins cannot diffuse on the surface. Due to the limitations of computational resources, only a small number of simulations were performed for with diffusion, but these confirm the fact that the corona has effectively reached its steady state within 360 s. Formally, we expect the rate-equation model in the limit of non-desorbing proteins to produce inaccurate results after an extreme period of time, as it can straightforwardly be seen that a steady-state solution with , = 0 can only exist when the surface of the NP is completely covered in proteins. Geometrically, this optimum packing cannot be achieved for spherical proteins. In practice, however, we find that for the set of proteins considered here the numerical integration of the rate-equation treatment provides physically realistic results even for timespans of up to one year, far exceeding experimental timescales.
For small but non-zero desorption rates, we find again that the evolution of the corona is significantly delayed. Even a relatively minor decrease in the desorption rates by a factor of 10 ( Figure 8) results in the steady-state corona taking several days to form in the hard-sphere model. A slightly worse agreement between the KMC and rate-equation models can be observed in the initial stages of Fig. 8. Taking into account the results found for irreversible diffusion, it is likely that this is due to the inability of the layer of bound proteins to restructure sufficiently quickly during the initial adsorption to provide the fluid-like properties necessary for the model to hold. We have performed further simulations of the initial stage of formation including diffusion and find a better agreement when this is included. The adsorption rate constants listed in Table 1 are quite low, corresponding to (on average) less than one potential adsorption event per second for Fib interacting with an NP of radius 10 nm, and even the highly abundant HSA undergoes less than one hundred events per second. To demonstrate the effect of either significantly increasing the concentration or mobility of proteins, we next investigate the outcome of increasing all adsorption rate constants by a factor of 1000. Although this produces unrealistically high collision rates, it serves to demonstrate an important point. The results of numerically integrating the rate equations for both the HS and MF models are shown in Fig. 9. In the MF model, there is essentially no effect on the evolution of the corona compared to the results obtained using the original rate constants. In contrast, however, we observe a significant difference in the HS model: the steady state corona is now primarily composed of the more abundant HSA rather than the more strongly binding HDL. Furthermore, the corona still takes a significantly long time to reach this steady state, with equilibrium not reached after 12 hours, although only minor changes take place after this point.
Due to the high rate of adsorption, simulating the corona using the KMC method for the timespans required to reach steady-state is impractical for this system even with diffusion disabled.   Adsorption rates for all proteins are increased by a factor of 1000 compared to the values presented in Table 1.
As a further test of the model, we add two additional proteins to the set originally considered and re-calculate the expected corona composition. The first of these two proteins, FP1, is simply the original set of parameters for Fib with the concentration increased to be equal to that of HSA, while the second additional protein FP2 is equivalent to HDL but with a slightly greater radius = 5.5 nm in comparison to the original radius of 5 nm. Results are plotted in Fig. 10 for both HS and MF models. In the HS model, it can be seen that even significantly increasing the concentration of Fib is not sufficient to cause it to be expressed in the long-term corona, as although it is initially adsorbed it is rapidly replaced by more strongly adsorbing proteins (HDL, FP2) or those with similar adsorption characteristics but a smaller size (HSA). Indeed, although both HSA and FP1 have nearly identical values of [ ] and , it is notable that of the two only HSA is exhibited in the steady-state corona. Meanwhile, FP2 has adsorption kinetics quite similar to that of the original protein HDL, but due to its increased size is less strongly expressed in the corona. Conversely, in the mean-field model, HDL and FP2 have identical surface coverages as a result of the fact in this model these coverages are independent of the radius of the protein. Thus, the difference in surface coverage observed in the HS model is due to the radius-dependent insertion probabilities. Likewise, in the MF treatment HSA and FP1 exhibit very similar profiles, in contrast to the size effects observed in the HS model.

DISCUSSION
From the results shown in the previous section, it is clear that the timescale for the formation of the corona varies significantly with the model used to describe how the rate of binding is reduced by pre-existing proteins. At one extreme, if proteins are allowed to move freely around the surface and deform to an effectively arbitrary degree as in the mean field model, the corona approaches equilibrium rapidly. At the other extreme in which the proteins are rigid and fixed in place, the time taken for the corona to reach equilibrium is dramatically increased for the set of proteins considered here. In both cases, the results of the KMC simulations can be adequately explained using analytical models of the adsorption process. Strictly speaking, the model of Ref. (27) used for the HS approach is valid for an equilibrium state, but appears sufficiently accurate for the dynamical states observed here. Qualitatively, both HS and MF models predict reasonably similar behaviour for the evolution of the corona -a sharp initial rise in HSA followed by slow replacement by HDL until the steady-state is reached, although on significantly different timescales. With the set of rate constants used here, the NP corona in the HS model takes several hours to reach equilibrium as a consequence of the decreased rate of adsorption at higher coverages. The difference in timescale from the HS approach arises purely from the decreased probability for proteins to find sufficient space on the NP to adsorb, and is independent of the rate constants. Indeed, even increasing the adsorption rate by a factor of one thousand still results in the HS model predicting timescales for corona formation on the order of twelve hours.
In contrast to the MF model, the HS model predicts that the size of the NP relative to the size of the proteins significantly alters the evolution of the corona. This is especially apparent in the results shown in Fig. 10, where simply increasing the radius of a protein by 0.5 nm significantly reduces the extent to which it features in the corona. Furthermore, it can be seen that in the MF model the more strongly-binding protein essentially entirely excludes a weaker-binding one, whereas in the HS model a weakly-binding but small protein may co-exist in the steady-state corona with a larger, more strongly binding protein. The question therefore remains as to which model -MF or HS -is more appropriate to describe the formation of protein coronas under biological conditions. In reality, proteins are neither completely deformable as assumed for the mean-field approach, nor completely rigid as in the hard-sphere model, and so we may expect that the actual timescales are not as slow as predicted in the HS model.
A key limitation of the analytical approach -whether mean-field or hard-sphere -is the fact that obtaining the solution to the set of differential equations rapidly becomes difficult as the number of proteins considered increases. For the results shown here, we have included a set of only three to five proteins, for which the equations may be numerically integrated using e.g. Mathematica. For the mean-field model at fixed concentrations of unbound proteins as considered here, an analytical time-dependent solution can be straightforwardly obtained using matrix methods, producing solutions which are sums of exponential functions. However, both the numerical and analytical approaches to solving these equations rapidly become more complex as the number of proteins grows and if we must also consider multiple orientations for each protein then even this simplified model may stretch the limit of computational resources, especially if the numerical integration must cover a large timespan. The KMC simulations, in contrast, scale reasonably well with respect to the number of proteins but become increasingly slow for larger NPs when a significant number of proteins are bound, especially if surface diffusion is enabled.
We must also address the limitations of the current model. We have assumed for the HS model that the proteins can be represented by hard disks with diameters calculated from that of the corresponding hard sphere. It has previously been shown in the RSA model that a true hard-sphere model exhibits different dynamics as a result of the decreased ability of smaller proteins to block larger ones, and, conversely, the ability of smaller proteins to bind in gaps between the NP and larger proteins (31). However, that work suggests that for proteins a hard disk model is more reasonable, and it is unlikely that proteins make only a single point contact with the NP as the true hard sphere model would require. For the system of proteins considered here, the sizes are sufficiently close that there is not likely to be a significant difference caused by the differences between these two models. A further complication is that protein shapes may deviate significantly from spherical and deform to improve their binding ability, with alternate binding profiles shown to result in increased coverage (38). While not accounted for in the present model, we note that the surface function (Eq. 6) derived from scaled particle theory is in principle valid for an arbitrary convex binding area. It is therefore likely that more complex protein shapes can be modelled by finding the area and perimeter of the projection of the protein onto the surface of the NP, treating different orientations as being effectively different proteins. This would further allow for a differentiation between different binding profiles, e.g. head-on vs. side-on for ellipsoidal proteins. Here, we do not investigate this further, but the CoronaKMC tool provided in the UnitedAtom package supports proteins with multiple binding profiles by representing these as hard-spheres of varying size.
The concentration of unbound proteins in this work and the rate at which they diffuse to the NP surface is here assumed to be fixed. Previous work in the MF model has indicated that finite rates of diffusion to the surface of the NP influence the dynamics of the corona formation (39), which could be represented here by a decrease in the adsorption rate . Likewise, it is reasonable to assume that the dynamics will be altered if the concentration of NPs [ ] is sufficiently high that substantial protein depletion in the environment occurs, i.e., if [ ] is of the same order of magnitude as the concentration of the protein. This latter effect can be be accounted for in the rate equation model by replacing [ ] → [ ] (0) − [ ] ( ), and in the KMC simulations by dynamically updating the concentrations of proteins based on the number currently bound to the NP. A further possible extension would be to include effects such as three-body interactions in the KMC simulations to produce better agreement with experiment as in Ref (21). Our model also does not capture more complex interactions between the proteins and their environment as considered in the dynamic density functional theory of Angioletti et al (24,25), which indicates that for an accurate model of corona formation it may be necessary to compute more physically realistic models of the protein density surrounding an NP, especially for NPs surrounded by a permeable gel. In the present case, however, we operate under the assumption that the NP and proteins in solution are sufficiently mobile that no significant concentration gradients develop.
It is important to note that the HS model, in which the proteins are represented as hard spheres, is conceptually quite similar to mesoscopic molecular dynamics simulations of the formation of protein coronas as in (21,22,40,41). In future work, it would be interesting to compare the predictions of this model to such simulations. However, the timescales predicted for the system to reach equilibrium are significantly longer compared to the timescales for individual events. Here, we assumed that proteins collide with the NP at a typical rate on the order of 1/second, and equilibrium takes ≈ 20000 seconds to reach. To get more realistic estimates of the corona formation times, the collision rates from real protein and NP concentrations should be used. On the other hand, it may therefore be possible to simulate only the initial formation of the corona using such simulations and require the HS model to extrapolate to greater timespans. A similar approach has been employed in Ref (21), where the MF model has been fitted to the results of the simulation to extrapolate to greater time periods than can be achieved through MD simulation. Based on the present work, in which the MF and HS models predict both very different timescales for formation based on the same rate constants and concentrations and differing steady-state corona compositions, we conclude that extraction of rate constants by fitting to numerical data (obtained through experiment or simulation) must be done with care to ensure that the correct underlying model is chosen. That is, attempting to fit the MF model to results obtained in the HS model would result in significantly different rate constants to the underlying values, and vice versa. Thus, it is vital that experimental results are interpreted by applying the correct model to obtain physically meaningful rate constants.

CONCLUSION
Based on the random sequential adsorption and desorption paradigm, we have developed a rate equation model for the evolution of the nanoparticle-protein corona for proteins modelled as hard spheres binding to to spherical and cylindrical NPs. We have demonstrated that the geometry of the NP directly influences the composition of the corona through altering the efficiency with which proteins can pack and new proteins can adsorb to the NP, and that this results in geometry-dependent corona compositions even when the rates at which proteins adsorb or desorb from surface sites are kept fixed. We have found that these effects can be explained using scaled-particle theory for both spherical and cylindrical NPs as long as there is either desorption or diffusion of adsorbed proteins to ensure a fluid-like surface layer. Our results pave the way for an improved understanding of experimental results and the ab initio prediction of protein corona compositions.

AUTHOR CONTRIBUTIONS
IR and VL designed the research. IR carried out all simulations, analyzed the data. IR and VL wrote the article.