Modelling Bi-specific Antibodies in Aqueous Solution

This study presents theoretical results for physico-chemical properties of system of molecules modeling bi-specific antibodies, such as, dual-variable-domain monoclonal antibodies (DVD-Ig) and Fabs-In-Tandem Immunoglobulin (FIT-Ig). These molecules are representatives of the engineered proteins that combine the function and specificity of two monoclonal antibodies. Individual molecules are here depicted as an assembly of nine (or in case of the Fit-Ig eleven) hard spheres, organized to resemble the Y-shaped object. The effects of the increased size, asymmetry, and flexibility of individual molecules on measurable properties of such systems of molecules are investigated. We examined the liquid-liquid phase separation, the second virial coefficient 𝐵 2 , and viscosity under various experimental conditions. The calculations are compared with the data for regular monoclonal antibodies and discussed in view of the experimental results for DVD-Ig solutions available in literature.


Introduction
Monoclonal antibodies (mAbs) are proteins of the immune system, engineered to suit therapeutic purposes [1][2][3] .They offer promising treatments for auto-immune disorders, cancer, and many other often fatal diseases.In developing of therapeutic proteins and their mixtures 4 the goal is to obtain formulations that are sufficiently concentrated in antibodies and yet free of protein aggregates, see Ref. 5 .The latter, not only limit their effectiveness as therapeutics, but may also pose a health risk for patients.
Single-target immunotherapy does not seem to destroy sick cells sufficiently.Binding multiple epitopes with a single (bispecific) antibody offers significant benefits, the most important one is much higher specificity.For review on possible future developments and application of these proteins see references [6][7][8][9][10] .In simple words we may say that bi-specific antibodies (bsAbs) combine the functionality of two antibodies in one molecule.
In literature we can find several Ig molecular platforms 7 , one is DVD-Ig and the other is Fit-Ig.This new generation of therapeutics can target two or even more disease mechanisms.They are designed for cancer immunotherapy 11 and currently several of them are in clinical use as also in diagnostics.They are designed as an Ig-like molecules, except that each light chain and heavy chain contains two variable domains.As such they are larger than an original immunoglobulins, asymmetric and in addition also more flexible 12 .Unfortunately, this increases their propensity for aggregation 3,10 .
In the present study we propose to analyse some measurable properties of the model DVD-Ig (dual variable domain immunoglobulin) [12][13][14][15][16][17] and FIT-Ig proteins in aqueous environment.These large antibodies (M  around 200 kDa and more), having increased asymmetry and flexibility, should be costly to simulate in detail.Actually, we do not know any such attempt.To examine the effects of the increased size and flexibility on thermodynamic properties we utilized here an improved version of the thermodynamic perturbation theory (TPT) for associating fluids 18 .
Original version of the theory [19][20][21] has been successfully ap-plied in several recent studies [22][23][24][25] .The theoretical improvements presented here, enables us to account for the difference in the bonding abilities of the square-well sites located on the surface of different hard-sphere monomers of the molecule.This difference appears as a result of blocking effects due to the rest of the hard-sphere monomers forming the molecule 18 .
As main results of our calculations we present the viscosity data for model proteins, the second osmotic virial coefficient and the liquid-liquid phase separation diagrams.We discuss the results in view of experimental data available in literature 15 and compare our findings with the results for regular mAbs studied recently [22][23][24][25] .

Models of the bispecific antibodies solution
The 9-and 11-bead representations of the antibodies studied here are modifications of the 7-bead model proposed and applied before [23][24][25] : the models are schematically presented in Figure 1.Such asymmetric (Y-like shape) "proteins" are constructed from nine or eleven (not seven as are regular mAbs) beads.For additional illustrations see also Fig. 1 of Ref. 7 and Fig. 1 of Ref. 12 .
We consider a solution of protein molecules, schematically shown in Figure 1, with the number density   .Molecules modeling antibodies are depicted by   tangentially bonded hard-sphere beads (see also Ref. 24 ), forming three-arm completely flexible star shaped molecules with   -bead arms attached to the central sphere.There is   =2 beads forming the base and   =3 and   = 4 beads for the two upper arms of DVD-Ig and FIT-Ig antibody models, respectively (Figure 1).Hard-sphere diameters of the beads composing the model protein are denoted   .The last two monomers of the three-and four-bead arms, which represent double variable domains, each bear one off-center square-well site, randomly placed on the surface.These sites are denoted as " and  ′ (for double variable Fab domain) and " and  ′ (for double variable Fab' domain).As before, there is only one off-center square-well site located on the terminal bead of the arm representing fragment crystallizable region and denoted as .The pair potential    (12) between these beads is where ℎ () is the hard-sphere potential, ,  = ",  ′ , ",  ′ , ,  and    (> 0) are the width and depth of the off-center squarewell sites, respectively, and  (12) is the distance between these square-well sites.Effects of the size, asymmetry and flexibility of such molecules on the physico-chemical properties of their solutions are investigated next.More precisely, we studied the liquid-liquid phase separation, the second virial coefficient,  2 , and viscosity, all under various experimental conditions.The results are compared with the data obtained for regular (7-bead) monoclonal antibodies and discussed in view of the experimental results for DVD-Ig solutions available in literature [15][16][17] .

Calculation of measurable quantities
Theoretical approach used here [19][20][21] has previously been explained in several papers, see, for example, [23][24][25][26] .Similarly as before, the Helmholtz free energy of the system  is presented as a sum of two terms, where     is the Helmholtz free energy of the reference system and Δ  the corresponding contribution due to protein association.Here the reference system is represented by the   -bead model with    = 0.According to our previous studies [23][24][25][26] we have where ℎ ( +  ) is the Percus-Yevick contact value of the hard sphere radial distribution function (RDF) and Δ ℎ is excess hard-sphere Helmholtz free energy.Δ  is calculated us-ing expression similar to that used earlier i.e.
where  takes the values  1 ,  2 ,  1 ,  2 ,  and   is fraction of the particles  bonded on a site .These fractions follow from the solution of the set of five equations where and (ℎ) ( +  ) is the contact value of site-site RDF between the centers of the hard-sphere beads bearing sticky sites of type  and  18,27-29 , i.e.
From here on, the chemical potential  and pressure  can be obtained via the standard relations, i.e.
and used for the phase equilibrium calculations.The densities of the coexisting phases follow from the solution of the set of two equations (, representing phase equilibrium conditions.Here,  ()  and ()  are the densities of low-density and high density phases, respectively.
Our calculation of the viscosity is carried out following the scheme developed earlier.According to this scheme 23 the relative viscosity of a solution, / 0 , can be given as where  is viscosity of the solution and  0 viscosity of the solvent, while  is the mass concentration of solute molecules (mg of pro-tein per mL of the solution).Further, (, ) is the weight fraction distribution of the clusters containing  molecules, while  () describes contribution of the cluster of the size  to the viscosity.For  () we use here the function applied previously 23 , i.e.
where  and  are the fitting parameters, obtained from the comparison of the theoretical and experimental data.We assume that the sites ", " and , which are located on the tips the variable and crystallizable domains, interact in-between with the same strength.Also, the sites  ′ ,  ′ , which are located on the monomers next to the terminal monomers of the variable domains, interact in-between and with the rest of the sites with reduced strength, i.e.    =  for ,  = ", ", ,    =   =  for  =  ′ ,  ′ and  = ", ",  ′ ,  ′ , .Notice that 0 ≤  ≤ 1.
The distribution function (, ) can be calculated using the information about fractions of the molecules   not boded through the site .For the model at hand we have a set of five fractions, which completely define (, ), i.e.   ′ ,   ′ ,   ′′ ,   ′′ and   .Note that due to the symmetry of the model we have: " =  " =   and   ′′ =   ′′ .Next we assume that distribution function (, ) can be approximated by the distribution function of the version of the model with three effective equivalent sites, i.e.

𝑃(𝑛, 𝛾)
where 30 and  is the fraction of the particles, which are not bonded via one of these sites.To relate both versions of the model we assume that this fraction is obtained from the equality of Helmholtz free energy for the current model and the model with three equivalent sites, i.e.
where for Δ  we use the expression (5) and Due to this relation thermodynamics of the model with equivalent sites is the same as that of the original model.Note that for =0 the five-site model reduces to the model with three equivalent sites.

Results and comparison with experimental data
In important experimental papers Raut and Kalonia [15][16][17] examined the liquid-liquid phase separation in solutions of DVDantibodies as also their viscosities.They measured the cloudpoint temperature,   , which is marking the onset of the phase separation.They have investigated the effects of size of the protein, the salt concentration, H, and nature of some other additives on the phase behavior and viscosity of the solution.Un- fortunately, our present model does not include electrostatics so it cannot account for these effects.We can therefor only investigate the influence of the increased size and asymmetry of the DVD molecule in comparison with regular antibody molecule studied before 23,24 .

Viscosity
We fit the model calculations to experimental data for the viscosity of antibody solutions presented in 16 .In that paper viscosity measurements for mAb and DVD-Ig protein solutions at different values of H, ionic strength  and protein mass concentration  =  2   /  were carried out.the molecular weights of mAb and DVD protein molecules were  2 = 150, 000 g/mol and  2 = 200, 000 g/mol, respectively, and   denoted Avogadro's number.The temperature of the solutions was 25 • C. The model parameters, which have to be fitted are: , , , ,  and .The value of  was as before chosen to be approximately equal to the hydrogen bond length, i.e.  = 0.15 nm.Further we assume that the size of hard-sphere monomers  for both 7bead and 9-bead models is the same and we are using here the value  = 2.5 nm.Also, in all the cases studied here we use the same value for the parameter  defining the strength of the cross interaction; =0.6.Parameters  and  describe contribution of clustering to viscosity of the solution.In this analysis we assume that these parameters depend only on the type of the molecules (either mAb or DVD-Ig) and are independent of the ionic strength of the solution and its H value.This is a severe approximation.
The choice of parameters used in calculations is collected in Table 1 and comparison of the experimental and theoretical results for viscosity are presented in figures 2-3.We note in passing that the chosen set of parameters is not unique.We have carefully examined the effects of reasonable variations of these parameters; small quantitative but not qualitative changes were noticed.
In Figure 3 we show our results and corresponding experimental data 17  conclusion is that the model can describe viscosity measurements reasonably well.Additional information comes from inspection of Table 1.Though our model does not account for electrostatics we can still learn some useful information from the variation of parameters in Table 1.It is interesting that  and  do not change with H of the solution neither with the added salt concentration so the results in figures 2 and 3 are fully determined by the parameter , i.e. by the protein-protein attraction.In other words, viscosity of the DVD-Ig solutions measured in Ref. 17 can be modelled with reasonable accuracy fixing values of all other parameters but the energy of patch-patch attraction, .Note that this parameter is defined as negative quantity (see Eq. 2) so smaller in magnitude value of  means weaker attraction between the protein sites (patches).For all the situations analyzed here the attraction decreases with an increasing salt content.This holds true for mAbs and DVD proteins.

Second virial coefficient
The second virial coefficient,  2 , quantifying the binary solutesolute interaction in dilute solutions, is one of the most important measurable quantities in protein solutions.It is known that value of this parameter can be used as an indicator of the crystallization 31 as also, that low  2 values are indicative for high viscosity 23,32 .This coefficient is defined as and can be obtained from the osmotic pressure equation as explained elsewhere 33 .It is most often presented as a function of the protein concentration   .In contrast to this here in Figure 4 we compare our predictions for the second virial coefficient of the DVD-Ig solution as a function of H, but at two different values of the ionic strength, i.e.
=15 and =50 17 .We found a qualitative agreement between our calculations and experimental data.At H=5.1 experimental measurements give small positive and small negative values for  2 at  = 15 mM and  = 50 mM, respectively.At the same time for this value of H and both values of the solution ionic strength theory predict almost the same small positive values, which are intermediate between those obtained experimentally.As H increases, both experimentally and theoretically calculated  2 decrease until around H≈6.1 for  = 50 mM and H≈6.5 for  = 15 mM, where  2 is negative and reach its minimum values.Thus at these values of H interaction between the protein molecules is strongly attractive.We need to know that p of the protein is around 7.5.This is also reflected in the behavior of the viscosity as a function of H, shown in Figures 2 and 3.Here the most rapid increase of the viscosity is observed for H=6.5.On the other hand while experiment predict almost the same values of  2 for solutions with  = 15 and  = 50 at H=6.1 (see Figure 4) corresponding values of the viscosity under the same conditions are quite different (see Figure 2 and Figure 6 of 17 ).Here viscosity of the solution with =15 is about two times larger than that with =50.This is reflected in the behavior of the theoretically calcu-lated  2 , i.e. here  2 ( = 15) ≈ 2 2 ( = 50).Further increase of H causes slight increase of the theoretically calculated  2 , which is still has a negative value.Similar behavior can be observed for the experimentally obtained  2 for the solution with =15 mM.For solution with =50 mM experimental  2 remains constant for H values in the range 6, 1 ≤ H ≤ 7.1 Thus the model used here enables us to reproduce in general correlation between viscosity and second virial coefficient of the system.For more accurate description of this correlation less coarse grained model will be needed.

Liquid-liquid phase separation
Here we investigate the effects of the size (molecular mass) and asymmetry on the liquid-liquid phase separation.We examined three different models presented in Figure 1(a), (b), and (c); i.e. the regular mAbs, and two bi-specific variants of antibodies (b) DVD-Ig as also and (c) FIT-Ig.
The calculations (see Figures 5 and 6) present data for the liquid-liquid phase separation modeling molecules shown in Figure 1.In this calculation all the proteins are interacting via the outermost beads: in case of (a)    =  (,  = , , ).For DVD-Ig and FIT-Ig model molecules (models (c) and (d)) we have:    =  (,  = ", ", ) and    =   =   = 0 (,  =  ′ ,  ′ ).Numerical results for this case are presented in Figure 5.In this figure we show our results using regular TPT1 approach 23 as also its improved version, the so-called modified TPT (mTPT) 18 , as discussed above.Note that within the mTPT approach, instead of the contact values of the hard-sphere RDFs, the contact values of site-site RDFs ( 9)- (11) were used.Here mTPT gives slightly wider phase diagrams with slightly smaller values of the critical temperatures.
In Figure 6 we present our results for the phase diagram of 9-bead DVD-Ig model with different values of the patch-patch interaction (measured by parameter ), which include the patches of the type  ′ and  ′ , i.e.
=  (,  = ", ", ) and Thus for =1 we have the 5-patch model and for  = 0 we have the 3-patch model.For the parameter  being smaller then one the sites are not completely coupled and the phase diagram gradually changes its shape.It becomes narrower around the critical temperature and wider for low temperatures, as it is shown in Figure 6.For ≈ 0.6 or less one can identify two distinct coexistence regions: at higher temperatures the phase diagram is narrow and coincide with the phase diagram for the three-patch version of the model and at lower temperatures the coexistence region is about three times wider and coincide with the coexistence region of the five-patch version of the model.Transition between these two regions occurs in a narrow window of the temperature with the width, which becomes smaller for smaller values of parameter .   =  ( ,  = ", ",  ) and   ′  =   ′  =   ( = ", ",  ′ ,  ′ ,  ) in  *   coordinate frame.Here from the top to the bottom at  = 0.06: =1 (green line),  =0.9 (blue line), =0.8 (pink line), =0.6 (red line) and  =0.4 (brown line).Black line represent results for =0.
phase behavior of the DVD-Ig protein solution at different values of H and ionic strength.The calculations were carried out for the same values of model parameters ,  and  as determined in our viscosity calculations above.The values for the potential well depth  are given in Table 1.However, the value of the potential well width  needed to reproduce the experimental results for the liquid-liquid phase diagrams of the solution at H=6.5 and two values of the ionic strength, i.e. =15 mM and =50 mM 15 , is in this case twice larger, i.e. =0.3 nm (see Figure 7).Unfortunately, the authors 15 did not provide the liquid-liquid phase separation graph in the form T* vs protein concentration be compared with calculations in the broader range.Their experimental data 15 are presented in Figure 4. Experiment predicted that there is no phase separation above  = 295 • K for the solution with =50 mM and above  = 298 • K for the solution with =15 mM, so we may take this value as an estimate of the crit- ical temperature.They also estimated for the critical density to be somewhere between 20 and 100 mg/ml.Theoretically calculated phase diagrams appear to be in a reasonable agreement with these values of the critical temperature and density.For the solutions with =50 mM and =15 mM critical temperature is   = 303.0• K and  = 318.4• K, respectively.In both cases the critical densities are the same, i.e.   = 30.8mg/mL.Similar values of the parameters were used to predict the phase behavior of DVD-Ig solutions at other conditions: H=7.0(=15 mM and =50 mM) and H=5.1 (=15 mM and =50 mM).These results are shown in Figure 8.Here location of the critical temperature for the phase envelopes at different values of H and  follows the general trend observed for the second virial coefficient (see Figure 4), i.e.  It has been shown by Bianchi and coworkers 30 that for patchy colloids the width of the LLPS envelope critically depends on the number of patches (i.e. the number of the off-center square-well sites) on the particles.In our one-component model the minimum number of sites fully bonded is three for the system to be able to form a network and to phase separate.However, by studying binary mixtures with varying the average number of attractive square-well sites 30 it has been possible to reduce the LLPS envelope width to very small values.At the same time the critical density was approaching zero.A similar effect has been later noticed studying the model mAbs molecules, confined in the hard-sphere fluid 34 .In our calculation above, such a non-complete binding is modelled by variable , varying the strength of the interaction.Experimentally, there are several ways to achieve this.The external parameters, which may be modified are: temperature, H, nature of the buffer and/or added electrolyte, as also the presence of other binding species (proteins).

Conclusions
In this paper we present a theoretical study of the model bispecific antibodies forming a liquid.The data published in Refs. 15,16provide useful guidance to behavior of the DVD Ig solutions, in particular with respect to the salt concentration, H, and nature of additives.Very valuable are data on the protein-protein interaction contained in the viscosity, the second virial coefficient, and the   measurements.
As already mentioned above our current model is not designed to capture the influence of H, concentration and the nature of the added buffer, as also other subtle effects.Despite of these shortcomings we have been able to qualitatively reproduce some important solution properties.In addition, we proposed an explanation for the very narrow width of the experimental liquid-liquid phase transition envelope.
As a weakness of our approach it might be considered the fact that we can model the viscosity behaviour with the particles having five fully bonding, while for the equally good agreement with the liquid-liquid phase separation data we need to make this interaction weaker ( smaller than 1).At this point we shall stress the facts that viscosity is a dynamic while LLPS is a thermodynamic property.In particular, viscosity depends on the state of solvent, so it is strongly temperature dependent.The viscosity measurement were taken at 25 • C, that is above the critical temperature.Taking into addition into account also the other factors which may influence the liquid-liquid phase separation (for example, the solution composition) the weakness mentioned above, is not that surprising.

Acknowledgments
TH and YVK acknowledge financial support from the National Research Foundation of Ukraine (Project no.2020.02/0317)

Fig. 2 Fig. 3
Fig. 2 Relative viscosity / 0 of mAb (dashed lines and empty triangles) and DVD-Ig (solid lines and filled triangles and circles) protein solutions as a function of protein mass concentration  at ionic strength  = 0 (upward triangles and black lines), =15 mM (downward triangles and red lines) and =50 mM (circles and blue line) at  = 25 • C and H=6.1.Symbols represent experimental 17 and and lines the theoretical (this work) results.

Fig. 4
Fig. 4 Second virial coefficient  2 vs H for the DVD-Ig solution at ionic strength  = 15 mM (black lines and symbols) and =50 mM (red lines and symbols).Here solid lines and filled symbols represent theoretical results and dashed lines and open symbols stand for the experimental results 17 .

Figures 7 and 8 Fig. 5
Figures7 and 8present our predictions for the liquid-liquid

Fig. 7 Fig. 8
Fig. 7 Liquid-liquid phase diagram for DVD-Ig protein solution in    coordinate frame at H=6.5 and ionic strength =15 mM (solid curve and filled circles) and =50 mM (dashed curve and empty circles).Here experimental results 15 are shown by the circles and theoretical results are depicted by the lines.Theoretical predictions for the critical temperature and density are shown by empty (=50 mM) and filled (=15 mM) triangles, respectively.

Table 1
Parameters of the models for mAb and DVD-Ig protein solutions