*2.2. Entropy Calculation*

For each pair of subsequent mers of GAG's chain, the time series (containing 1000 points) of the dihedral angles have been obtained. As described in Section 2.1, following pairs of these angles: Φ1−<sup>4</sup> vs. Ψ1−<sup>4</sup> and Φ1−<sup>3</sup> vs. Ψ1−<sup>3</sup> have been analyzed [31]. 2D histograms have been calculated from angles values of all subsequent pairs of mers (and all time steps). Histograms of the stronger bound structures (from the 10 picked as described in Section 2), for different ionic solutions, have been presented in Figure 3 for CS6 and in Figure 4 for HA.

Following an approach described in [44] (see Figure 8 therein), the conformation entropy has been computed from the 2D histograms of pairs of angles. In more detail, from each histogram, the Shannon entropy [22,23,30] has been computed using the formula:

$$S = -R\_0 \sum\_{i,j} p\_{i,j} \log(p\_{i,j}).\tag{2}$$

Here, *pi*,*<sup>j</sup>* is the empirical probability of the first angle being in the *i*th bin and the second angle being in the *<sup>j</sup>*th bin, and *<sup>R</sup>*<sup>0</sup> = 8.314 J·K−1·mol−<sup>1</sup> is the (scaling) gas constant.

**Figure 3.** Normalized histograms of values of angles (Φ1−4, Ψ1−4) (top) and (Φ1−3,Ψ1−3) (bottom) for main chain CS6 with (**a**) Na+, (**b**) Ca2+, and (**c**) Mg2+. Angles were taken from the whole time series of the YASARA simulation. The symbol *n* is a number of angles' pairs and is equal to number of angles type (24 for angles 1, 4 and 23 for angles 1, 3, cf. Section 2.1) multiplied by number of time points (1, 000).

**Figure 4.** Normalized histograms of values of angles (Φ1−4, Ψ1−4) (top) and (Φ1−3,Ψ1−3) (bottom) for main chain HA with (**a**) Na+, (**b**) Ca2<sup>+</sup>, and (**c**) Mg2<sup>+</sup>. Angles were taken from the whole time series of the YASARA simulation. The symbol *n* is a number of angles' pairs and is equal to number of angles type (24 for angles 1, 4 and 23 for angles 1, 3, cf. Section 2.1) multiplied by number of time points (1, 000).

## **3. Results**

Some HSA segments are more prone to creating intermolecular interactions than others. The complexity of protein–GAG interactions is in part caused by the conformational flexibility of the GAG's chain. An affinity of the HSA to GAGs, firstly tested by docking method, has been present in Table 1 sorted by binding energy. While the docking method relies on adjusting a ligand to a receptor in crystal form, then putting the complex into a water solution changes the intermolecular interactions map. After equilibration and 100 ns of MD simulations, the binding energy changed, and the order of best-bound complexes changed. In the case of CS6, the new order depended on added ions, and its value (averaged over three realizations' binding energies with a different salt added) has been written in the first column of Table 1 in the brackets. HSA binding sites did not change much during the MD simulations. However, the number of interactions such as hydrophobic-polar, hydrogen bonds, ionic, and bridges have changed.

**Table 1.** Binding ranks of HSA-CS6 complexes. The first column contains: rank after docking (averaged rank after MD simulations). The strongest connected domains are marked in bold letters.


A closer analysis of the interactions for HSA-CS6 complexes will be the subject of another study, similar to the ones performed for HSA-HA complexes [15]. In the present paper, we are focused only on the conformational entropy of the GAGs chains. In the case of CS6, the most stable turned out realization was number 2, thus the situation when CS6 wrapped around the HSA and bound to IA-IB-IIA-IIIA-IIIB domains had the strongest binding to IIIA. No matter what ions have been added to the system, the binding energy stays high compared to the rest of the realizations (thus with different initial conditions of the binding map). After docking, amino acids that created the higher number of interactions with CS6 were Glu and Thr, and next in frequency: Lys and Asp. MD simulations show that Arg and Lys are more prone to create more ionic interactions and hydrogen bonds than the other amino acids due to their positive charge with negatively charged sulfate and carboxyl groups.

In the case of HA, a similar situation can be seen. The most stable was complex 2, i.e., composed of HSA's domains IA-IB-IIIA-IIIB. The most binding amino acids, in general, were Thr, Glu and Lys. It is very important that domains IB, IIIA, and IIIB are key domains for the HSA transport function responsible for the heme-binding site (IB), Sudlow's site II (IIIA), and thyroxine-binding site (IIIB) [45]. Comparing Table 1 to Table 2, one can notice that the HSA best binding segments differ slightly between HA and CS6 but are similar for the two first best-docked structures. In both cases, IIIA and IIIB domains clearly prevailed in creating the highest number of interactions between GAG and the protein. This is because domains IIIA and IIIB are domains that show a positive net charge on the surface that allows for binding with negatively charged GAGs. In addition, in both cases, a fragment of GAG's chain strongly interacted with an IA domain. Analyzing differences between the values of energy of binding (averaged) in both cases, HA bound to albumin about 10% stronger than CS6 (see Table S1 in Supplementary Materials).

**HSA-HA Complex Rank HSA Binding Sites** 1(4) IA-IB-IIIA-**IIIB** 2(1) **IA**-IB-IIIA-IIIB 3(6) **IA**-IB-IIIA-IIIB 4(10) IIIA-**IIIB** 5(5) IIB-**IIIA**-IIIB 6(8) IA-IIIA-**IIIB** 7(2) **IA**-IB-IIIA-IIIB 8(9) IIIA-**IIIB** 9(7) **IA**-IB-IIIA-IIIB 10(3) **IIA**-IIB-IIIA

**Table 2.** Binding ranks of HSA-HA complexes. The first column contains: rank after docking (averaged rank after MD simulations). The strongest connections are marked in bold letters.

The method for computation of conformational entropy, based on a Ramachandrantype plot created for the pairs of dihedral angles (Φ,Ψ) [31,42], has been discussed in Section 2. Results of the computation have been presented in Figure 3 for the CS6 and in Figure 4 for the HA. Both of the results have been presented only for the best-bound complexes because the rest of the results had very similar characteristics. For comparison, three different realizations of YASARA simulations (thus, realizations with different initial structures) have been presented in Supplementary Materials Figures S1 and S2. The most probable angles, taken from Figures 3 and 4, have been presented in Table 3.

In the case of CS6, (Φ1−4, Ψ1−4) angles arranged in few clusters at ranges about: −150◦ – −60◦ for Φ1−<sup>4</sup> and −180◦ – −60◦ for Ψ1−<sup>4</sup> with the highest probability of occupancy near two more narrow angle ranges with a maximum at (−72◦, −76◦) (cf. Figure 3, top line). There can also be seen, however, the second angle region of about: 50◦–180◦ for Φ1−<sup>4</sup> and −180◦ – −50◦ for Ψ1−4, also with few narrowed clusters with a high probability of occupancy. The places of the spots and their intensity differ slightly between simulations with different ionic solution. In particular, in the case of Ca2+, most of the angles have been centered in one specific range around (−104◦, <sup>−</sup>76◦), while, in the cases of Mg2<sup>+</sup> and Na<sup>+</sup>, more than one high probability place can be seen. In the case of Mg2<sup>+</sup>, the distribution of the angles is the most uniform but with, similar to the Na<sup>+</sup> case, a maximum at (−68◦, −76◦).

**Table 3.** Most frequent (Φ1−4, Ψ1−4) and (Φ1−3, Ψ1−3) angles for the best bound complexes of HSA-GAG.


The Ca2<sup>+</sup> cations are distinguished by the fact that they form many more ionic interactions with the molecules than Na<sup>+</sup> and Mg2<sup>+</sup>. This can be the reason for the differences in the (Φ1−4, Ψ1−4) distribution plot. The greater probability cluster, which contains the maximum, is placed in the range of angles obtained lately in [46] for similar computer simulations of various kinds of chondroitin sulfate in water solution without any protein contribution. Thus, it can be clearly seen that the vicinity of HSA changed this crucial angle distribution, making them more disordered but still in a specific way. For proteins, similar (Φ, Ψ) angles are responsible for the formation of right-handed *α*-helices [42].

Despite the different chemical nature of protein and GAG molecules because peptide groups are linked at an *α*-carbon atom, not an oxygen like in the case of GAGs, the dihedral angles show how the chain-building units are rotated across the whole chain. Thus, the output secondary structure looks similar. In [46], there is, however, a lack of the second region on the angles, observed in our case with the positive values of Φ1−<sup>4</sup> angles.

In the case of angles (Φ1−3,Ψ1−3) for CS6, we can see quite different plots than the ones for (Φ1−4,Ψ1−4) (see Figure 3 bottom line, and Figure S1 in Supplementary Materials). According to [46], the differences were expected because, for angles 1–3, the most probable occupancy should be in the region of −100◦–−30◦ (Φ) and 70◦–180◦ (Ψ). In our results, the shadow (slightly visible red color on the plot) of those angles can be seen, especially in the case of Na<sup>+</sup> and Mg2<sup>+</sup>, where the maximum is placed within this range (cf. Table 3), near (−76◦, 169◦). Most of the angles have had the values similar like in the case of (Φ1−4,Ψ1−4) but more focused on regions near (−86◦, <sup>−</sup>76◦) (maximum for Ca2<sup>+</sup> case).

Conformational entropy for most cases is in the same range. However, there are noticeable differences in entropy between 1–3 and 1–4 angles. HA 1-4 angles show lower entropy i.e., are more stable than 1–3. The opposite is true for CS6. This can be explained by neighbouring groups. Carboxyl in HA and sulfur for CS6 form more stable contacts. Although acetyl group is highly reactive, it does not influence stability as dominant groups. The carboxyl group in CS6 is still weaker as compared to sulfur, which makes the contact more stable, as can be seen in Figure 5. The same behavior can be seen in Figure 6. The introduction of divalent ions increases entropy due to their destabilizing impact on protein, which is even more prominent for concentrations used in the present study.

Molecular conformational space available for HA chain in solution has been studied in [31]. The authors have searched for stable ordered forms of HA and have found many helices-type conformations (right- and left-hand side) that the HA chains prefer. Their findings based on potential energy computation for specific (Φ1−4,Ψ1−4) and (Φ1−3,Ψ1−3) angles are presented in Figure 2 of [31] . They obtain the location of three distinct regions with minimum potential energy surfaces. A main region consists of two wells (denoted as A–B in Figure 2 of [31] mentioned). Using the approach presented in this study, it is possible to compare the preferred dihedral angles for the HA chain, which is placed alone in the solution and in the vicinity of the HSA protein. In general, all these A-E regions (cf., Figure 2 in [31]) have been found in presented simulation results, but the intensity of these regions on the probability map varies depending on realizations (thus binding sites) and ion addition (see Figure 4 and Figure S2 in Supplementary Materials. In [31], for (Φ1−4,Ψ1−4) angles, the main A–B region is placed in the area about −120◦–−60◦ for Φ1−<sup>4</sup> and −180◦–−100◦ for Ψ1−4, which is in accordance with results obtained in Figure 4 for the best docked HSA-HA complex. The probability of finding the angles in clusters A and B is almost equal in cases of realizations with the addition of Na<sup>+</sup> and Mg2<sup>+</sup>, but, in the case of Ca2<sup>+</sup>, B prevails over A. The maximum has been found about (−72◦, <sup>−</sup>126◦) (B) and (−97◦, −155◦) (A). Moreover, a few different lighter clusters have been found similar to regions C, D and E.

In Figure 4, one can see two clearly identified clusters for (Φ1−3,Ψ1−3) angles and two smaller ones. The first of the bigger cluster, with a maximum at about (−115◦, −61◦) in the case of Na+, suits region C in [31], and the second, near the maximum (−54◦, 151◦) in the case of Mg2<sup>+</sup>, is placed in region A overflowing to region B. Region D, about (50◦, 120◦) angles, is also present.

**Figure 5.** Median, maximal and minimal entropy (over *N* = 10 realizations) for chosen ions and angles pairs for CS6 **(left panel)** and HA (**right panel**). As we have *N* = 10 realizations, the minimal entropy can by considered as the estimate of the 10'th percentile (first quantile), and the maximal one as the estimate of the 90'th percentile (9'th quantile).

Median, minimal and maximal values (over realizations) of entropy for various ions, angles and GAGs are presented in Figure 5. Bear in mind that, as we used *N* = 10 realizations, the minimal value can be used to roughly estimate first quantile, while the maximal value to roughly estimate the 9'th quantile.

In Figure 5, one can observe that, for analysed angles, ions and GAG type, entropy was in general greater for (Φ1−4,Ψ1−4) angles than (Φ1−3,Ψ1−3) in the case of CS6, but it was the opposite situation for HA: the entropy was slightly greater for (Φ1−3,Ψ1−3). As lower entropy shows lesser disorder in the system, the most stable systems were those with Na<sup>+</sup> ions added in both cases: HA and CS6 (cf., Figure 5). Complexes with Ca2<sup>+</sup> ions usually have had slightly higher entropy. Referring to the informative interpretation of entropy, one can conclude that, in the case of CS6, the pair (Φ1−3,Ψ1−3) carries significantly less information of the system than the pair (Φ1−4,Ψ1−4). This is not the case for HA. The difference is related to the presence of the sulfate group in GalNAc in CS6 that is more prone to forming hydrogen bonds and ionic contacts with HSA.

Entropy values for CS6 and HA with different ions, taken separately for each of computer experiment realizations, have been presented in Figure 6.

**Figure 6.** Conformational entropy for CS6 with (**a**) Na+, (**b**) Ca2+, and (**c**) Mg2+, and HA with (**d**) Na+, (**e**) Ca2<sup>+</sup> and (**f**) Mg2<sup>+</sup>.

Relatively large variations of the entropy between realizations are observed in the case of HSA-HA simulation results. Entropy varies within the range of 53–62 <sup>J</sup> Kmolfor CS6 and within the range of 55–61 <sup>J</sup> Kmol for HA. Thermal noise, or some non-equilibrium processes, can have some effect on these variations, but also the different binding sites of the two molecules (that varies between realizations) are of big importance for the entropy behaviour. The variations of entropy values may also coincide with rather high estimation error of this value.

The hypothesis that entropy value is tied to the value of the binding energy between the protein and the GAG is not supported by our simulation results. In more detail, the smallest entropy value, roughly 53 <sup>J</sup> Kmol , has been obtained for realization number 8 of the HSA-CS6 complex with Na+. This was the case, where only IA and IB domains of HSA were bound to CS6; thus, the protein did not affect the conformation of the GAG's chain by deformation of the (Φ,Ψ) angles much. The highest entropy value has been reported for case 3 of the HSA-CS6 complex with Mg2+. Here, the binding site was very similar to the one with greater entropy (number 8): IA and IIA. Entropy records for realization number 2 of the HSA+HA complex with addition of Ca2<sup>+</sup> are very interesting. In this realization, there were a huge difference between entropy for (Φ1−4,Ψ1−4) angles and (Φ1−3,Ψ1−3) angles. In this case, the HSA and HA molecules were best bound after MD simulation from all the realizations.
