Analysis of stabilization mechanisms in -lactoglobulin-based amorphous solid dispersions by experimental and computational approaches

Our previous work shows that β -lactoglobulin-stabilized amorphous solid dispersion (ASD) loaded with 70 % indomethacin remains stable for more than 12 months. The stability is probably due to hydrogen bond networks spread throughout the ASD, facilitated by the indomethacin which has both hydrogen donors and acceptors. To investigate the stabilization mechanisms further, here we tested five other drug molecules, including two without any hydrogen bond donors. A combination of experimental techniques (differential scanning calorimetry, X-ray power diffraction) and molecular dynamics simulations was used to find the maximum drug loadings for ASDs with furosemide, griseofulvin, ibuprofen, ketoconazole and rifaximin. This approach revealed the underlying stabilization factors and the capacity of computer simulations to predict ASD stability. We searched the ASD models for crystalline patterns, and analyzed diffusivity of the drug molecules and hydrogen bond formation. ASDs loaded with rifaximin and ketoconazole remained stable for at least 12 months, even at 90 % drug loading, whereas stable drug loadings for furosemide, griseofulvin and ibuprofen were at a maximum of 70, 50 and 40 %, respectively. Steric confinement and hydrogen bonding to the proteins were the most important stabilization mechanisms at low drug loadings ( ≤ 40 %). Inter-drug hydrogen bond networks (including those with induced donors), ionic interactions, and a high T g of the drug molecule were additional factors stabilizing the ASDs at drug loading greater than 40 %.


Introduction
An amorphous solid dispersion (ASD) is a molecular mixture of drugs and excipients that potentially increases dissolution rate and apparent solubility, and hence the bioavailability of poorly water-soluble drugs.(Vasconcelos et al., 2007;Schittny et al., 2020) However, the inherent propensity of the amorphous form to convert into the crystalline one poses a challenge because bioavailability relies on the specific type and concentration of the excipient.Recently, ASDs based on β-lactoglobulin (BLG) have been used for amorphous stabilization, dissolution, and solubility enhancement.For example, the incorporation of amorphous indomethacin in BLG keeps ASDs at drug loadings from 10 % to 70 % in an amorphous state at 40 • C under dry conditions for at least 12 months.(Kabedev et al., 2022) In another study, amorphous nimodipine stabilized by BLG exhibited better physical stability than with poly(vinylpyrrolidone) (PVP)-based ASDs at 10 and 17.5 % drug loadings.(Leng et al., 2023) Furthermore, their solubility was significantly higher than the respective nimodipine-PVP ASDs, in both simulated gastric fluid and simulated intestinal fluid.Finally, a rifaximin-BLG ASD showed a faster dissolution rate and higher drug solubility in three release media with different pH values (1.2, 4.5, and 6.5) compared to polymer-based ASDs and nanocrystalline rifaximin.(Zhuo et al., 2023) The BLG ASD approach has good efficacy, but the exact mechanism behind the stabilization of amorphous drugs in the solid state and the solubility improvement is still unclear.In an earlier study, the stabilizing mechanism of BLG for indomethacin in ASDs was investigated by combining differential scanning calorimetry (DSC) and molecular dynamics (MD) simulations.Samples formed single-phase ASDs with one glass transition temperature (T g ) at drug loadings below 40-50 %.For drug loadings above 40-50 %, the MD simulations revealed systems with a drug-rich phase.This finding was further supported by the identification of two T g events in the DSC measurement.Additionally, BLG-stabilization of indomethacin was at least as efficient at drug loadings of up to 70 %.
The stabilization mechanism was attributed to a combination of two factors: the reduced mobility of the indomethacin molecules in the first drug shell and the presence of hydrogen bond networks in the second and subsequent shells.The former factor dominates at lower drug loadings (up to 30-40 %), where indomethacin molecules are situated on the surface of BLG and surrounded by other BLG molecules, restricting their ability to recrystallize.At loadings greater than 40 %, hydrogen bond networks between indomethacin and BLG extend into the second and subsequent drug layers.This hydrogen bond network effectively stabilizes the entire mass at drug loadings of up to 70 %.(Kabedev et al., 2022) The indomethacin-BLG ASD study sheds some light on BLG-based ASD loading and stabilization mechanisms.However, only a limited number of drugs have been tested in the system, so it is not clear which drug properties are the crucial ones for successful preparation and stabilization of protein-based amorphous formulations.Among other things, the influence of melting temperature (T m ), pKa, T g , and hydrogen bond acceptors/donors has not been properly investigated.Using gelatin 50PS as a carrier to formulate solid dispersions, van den Mooter et al. studied twelve compounds, differing widely in their physicochemical characteristics.Their results revealed that eleven of these drugs could form true glass solutions at lower drug loadings.However, when the drug content was further increased, up to 40 % for some compounds, true glass solutions could not be achieved when processed by lyophilization (Pas et al., 2018).It is therefore essential to identify the potential effects of drug physicochemical properties on the formation and stabilizing mechanism of protein-based ASDs, in order to develop optimized formulations.
Prior findings indicate that the formation of hydrogen bonds is a crucial factor in stabilizing indomethacin in BLG-based ASDs.Indomethacin has one donor and four acceptors evenly distributed across it.As such, the prerequisites for the formation of a hydrogen bond network are met.However, this stabilizing mechanism can be absent, or at least not as efficient, for drugs with either only donors or acceptors, or for drugs with donors and acceptors if they are all located in the same part of the molecule.To this end, we spray-dryed BLG-based solid dispersions, loaded with five model drugs, each exhibiting different combinations of characteristics, including the number of possible hydrogen bonds.As shown in Table 1, these drugs also feature different acidity and alkalinity properties, as well as solid state characteristics as identified by the T m and T g values.Using a combination of experimental techniques and computer simulations, we investigated this panel of drugs to obtain a more comprehensive overview of the stability of BLG-based ASDs.In a series of computer simulations, we analyzed the diffusivity of the molecules, the propensity of the APIs to form hydrogen bond networks, and their affinity to the proteins.Finally, we developed an MD search strategy for the identification of crystalline patterns that might lead to recrystallization of the formulations.

Preparation of ASDs by spray drying
BLG-based ASDs samples were prepared with a Büchi B-290 spray dryer (Büchi Labortechnik AG, Flawil, Switzerland) equipped with an inert loop (Büchi B-295) and dehumidifier (Büchi B-296) in a closed configuration.All five drugs and the BLG were dissolved in 85 mL of ethanol and 10 mL of water, respectively.Then 5 mL of glacial acetic acid was added to the water and stirred at 300 rpm for 20 min.Subsequently, the BLG solvent was added to the ethanol containing the dissolved drugs to obtain a mixed solvent (final concentration 85 % ethanol; 10 % water and 5 % glacial acetic acid).The amount of dissolved solid content corresponded to a total of 625 mg of drug and BLG (from 10:90 to 90:10 in 10% increments for rifaximin, ibuprofen and ketoconazole; from 10:90 to 80:20 in 10 % increments for griseofulvin).Due to the limited solubility of griseofulvin in ethanol, a reduced total of 312 mg of (drug and BLG) was used at a 90:10 ratio (griseofulvin: BLG).
After stirring for 20 min, the solutions were spray dried with the following conditions: inlet temperature, 100 • C; outlet temperature 50-57 • C; feed rate, 10 mL/min; atomization air flow rate, 473 L/h: drying air flow rate, ca.35 m 3 /h.

X-Ray powder diffraction
Diffractograms of all samples were recorded using an X'Pert PANalytical PRO X-ray diffractometer (PANalytical, Almelo, The Netherlands) with CuKα radiation (λ = 1.54187Å).All samples were scanned in reflection mode from 5 • to 30 • 2θ, with a scan speed of 0.067 • 2θ/s and a step size of 0.026 • The acceleration voltage and current were 45 kV and 40 mA, respectively.The data were collected and analyzed with the software X'Pert Data Collector (PANalytical).XRPD diffractograms are presented in Fig. S1.

Volatile content determination
The volatile content of the freshly prepared formulations (moisture or residual solvents) was determined using a Discovery thermogravimetric analyzer (TGA) from TA Instruments Inc. (New Castle, DE, USA).Approximately 10 mg of powder was placed in an open platinum pan and analyzed at a heating rate of 10 • C/min up to 300 • C under 25 mL/ min nitrogen gas purge.The volatile content was determined as weight loss between 25 and 150 • C (n = 1).

Differential scanning calorimetry
The glass transition temperature (T g ) of the freshly prepared formulations was measured with the differential scanning calorimeter (DSC) from TA Instruments.All experiments were performed under 50 mL/min nitrogen gas purge.Approximately 6-8 mg of powder was weighed into aluminum Tzero pans and compacted using a cylinder to a : hydrogen bond donors.b : hydrogen bond acceptors;.c data from (Gushurst et al., 2014).d data from (Blaabjerg et al., 2017).e data from (Blaabjerg et al., 2017).f data from (Al-Obaidi and Buckton, 2009).g data from (Fung et al., 2018).
X. Zhuo et al. ensure homogeneity.Subsequently, samples were sealed with perforated hermetic lids except for the ibuprofen-BLG ASDs (sealed with intact hermetic lids).All samples were measured under the modulated temperature mode as follows.Samples with perforated hermetic lids were determined with a heat-cool-heat cycle, which was first annealed at 125 • C at a heating rate of 10 • C and then kept isothermal for 1 min to remove any residual moisture/solvent before being cooled to 20 • C. Subsequently, they were heated to 250 • C at a heating rate of 2 • C/min with a modulated temperature amplitude of 0.2 • C and a period of 40 s (n = 2).The samples with intact hermetic lids were kept isothermal at 0 • C for 10 min and then heated to the maximum of 250 • C at 2 • C/min with a modulated temperature amplitude of 0.2 • C and a period of 40 s (n = 2).The T g was defined as the midpoint over the step change of the reversing heat flow signal using the TRIOS software (version 4.1.1).
Representative DSC thermograms of freshly prepared ASDs are presented in Fig. S2.

Physical stability
All ASDs samples were stored in desiccators over silica gel at 40 • C under dry conditions.Each sample was regularly subjected to X-ray power diffraction (XRPD) analysis to detect any potential recrystallization.Stability assessment of the samples was done for up to 12 months (see Table S1).

Molecular topology parameterization
The generalized Amber force field (GAFF) was used in all-atom MD simulations.(Wang et al., 2004) The BLG model was developed in the course of the previous study on indomethacin ASDs with software packages Stage and Modeller.(Lundborg and Lindahl, 2015;Šali and Blundell, 1993) The topologies of griseofulvin and ketoconazole were developed earlier in our group and used in a study on the stability of HPMC-AS ASDs.(Edueng et al., 2021) Furosemide, ibuprofen, and rifaximin topologies were developed in the current study using the Stage software.

Simulation setup and parameters
ASDs were simulated with drug loadings of 40, 60 and 90 %.To satisfy the minimum-image convention at the lower drug loadings, eight BLG molecules were added to all systems.Each of the eight molecules was randomly rotated to avoid artifacts caused by artificial alignment of the proteins.Afterwards, API and water molecules were added to the box to form corresponding API:BLG proportions and to have 5 % water representing residual moisture (Kabedev et al., 2022).Packing of the molecules was done with Packmol (Martínez et al., 2009).Periodic boundary conditions were applied in three dimensions in all simulations.
As in the original study setup (Kabedev et al., 2022), the assembled systems first underwent 4500 steps of steepest descent energy minimization, followed by NVT (Nosé− Hoover thermostat, τ t = 2 ps) and NPT equilibration stages (Berendsen barostat, τ p = 5 ps, with a compressibility of 4.5 × 10 − 5 bar).(Nosé, 1984;Hoover, 1985) At the production stage, the v-rescale thermostat was used for temperature coupling, and the Parrinello− Rahman barostat for isotropic pressure coupling.(Parrinello and Rahman, 1981) After the equilibration stage, the systems were annealed.First, the temperature in the systems was raised from 298 to 510 (furosemide), 520 (ibuprofen), or 650 K (griseofulvin, ketoconazole, rifaximin) during 20 ns, held at that level for another 80 ns, then cooled to 200 K over the next 20 ns, before finally being raised back to 298 K within 20 ns.The total run-time for each simulation was 200 ns, and during the final 60 of these 200 ns, the temperature was maintained at room temperature.This thermal annealing was necessary to properly mix the systems and remove the memory of the initial random configuration.In the original study on indomethacin-BLG ASD, no significant morphological changes of the proteins were found during the annealing process, for which the maximal temperature was 550 K and held 60 ns.Here, we raised temperatures up to 650 K and held them at maximum for longer in order to melt drugs with higher T m (griseofulvin, ketoconazole, rifaximin).In order to ensure the original state of BLGs, we introduced restraints for the proteins in all ASDs during the first 200 ns.Thereafter, another 200 ns of simulation was added without restraints at room temperature, of which the final 150 ns were used for the analysis (see the timeline in Figure S3).This accelerated the conditions, BLG denaturation was avoided, and no artificial restraints were present under the final production run.The timestep during the production stage was set to 2 fs.Each system was run in duplicate with different, initial random placement of the molecules.The P-LINCS algorithm (Hess et al., 1997) was used to constrain bonds involving hydrogen atoms, and particle mesh Ewald (PME) summation (Darden et al., 1993) was used for electrostatic interactions.Cut-off distances of 1.2 nm were set for both electrostatic and Lennard-Jones interactions.Gromacs version 2018 was used for consistency with the previous set of simulations on indomethacin-BLG ASDs (Abraham et al., 2015).The TIP3P model was used as a water model (Jorgensen et al., 1983).

Analysis
We analyzed the data with Gromacs built-in tools, (Abraham et al., 2015) VMD software, (Humphrey et al., 1996) and in-house scripts written in bash, Tcl, Python, and MATLAB (version R2022b).The code of the most analytical scripts is available on GitHub and is thoroughly commented (github.com/computationalpharmaceutics).The measurements for all analysis routines were done at drug loadings of 90, 60 and 40 % in the equilibrated state and were averaged over simulation duplicates.For better readability we used the term ASD for all formulations studied computationally, including those that did not remain amorphous in the experimental part of the study.

Diffusivity of drug molecules.
High mobility of API molecules in an ASD is an indicator of low stability of the formulations.(Brunsteiner et al., 2018) It is often quantified in molecular dynamics with molecular diffusivity.We used the built-in Gromacs gmx msd tool to measure the diffusivity of each individual molecule within the last 150 ns of the production stage trajectories.Average values, standard deviations, and number of molecules were used to estimate the overall average and standard deviation of the API molecules for the specific drug loadings.Mean square displacement graphs of ten randomly chosen molecules in each API-DL combination were linear throughout the final 150 ns of the simulations.
To visualize the most mobile molecules, we used an in-house Tclscript to set the transparency of each API molecule in VMD, based on its diffusivity as measured with Gromacs.Two thresholds were set: 5 × 10 − 10 cm 2 /s and 5 × 10 − 9 cm 2 /s.At these values, drug molecules with diffusivity values lower than the former one would be completely transparent and those with diffusivity higher than the latter would be completely opaque.The range limits were chosen to ensure the informativity of the figures.Thus, the few mobile molecules present in all systems were omitted so that the few, highly mobile molecules would not skew the color scaling of the rest.Molecules in contact with the protein were also differentiated from those more distant from the BLGs to render the molecular mobility distribution more clearly.

Affinity and distance of APIs from the protein.
Affinity of the drugs for the protein contributes substantially to the stability of the ASDs.Any tendency of the drug molecules to self-aggregate in the formulation might limit the ASD efficiency, so it is important to take BLG-API interactions into account.While affinity alone is insufficient for predicting solid dispersion stability, (Yani et al., 2017) together with other measurables it partially clarifies the stabilization mechanisms.Here, we used a simplified approach for measuring the number of BLG-BLG, BLG-API, BLG-water, and BLG-ions contacts in different ASDs X. Zhuo et al. and presented them as a portion of the total contacts with all molecules.We could thereby estimate the tendency of the proteins to aggregate in ASDs with specific APIs.
We also measured the minimum distance between each API molecule and the BLG protein surface, as a function of time.Likewise, we measured the correlation coefficient between such minimum distances and diffusivity values of respective drug molecules.Both the minimum distance and number of contacts were measured with Gromacs tool gmx mindist.
2.3.3.3.Hydrogen bonds.One of the motivations for performing this study was to analyze the extent to which hydrogen bond networks contribute to BLG-ASD stabilization.We analyzed the average number of hydrogen bonds formed between both API-API molecules and between APIs and BLGs.A correlation analysis between molecular diffusivity and the presence of both types of hydrogen bonds was also performed for individual molecules.Hydrogen bond analysis was done with gmx hbond and in-house scripts.

Crystalline pattern search strategy. Pre-crystalline nuclei in
ASDs indicate potentially unstable formulations.We used information from the Cambridge Crystallographic Data center (CCDC) (Groom et al., 2016) to find crystalline patterns for furosemide (CSD entry: FUR-SEM01) (Lamotte et al., 1978), griseofulvin (CSD entry: GRISFL07) (Su et al., 2018), ibuprofen (CSD entry: IBPRAC03) (Stone et al., 2009), ketoconazole (CSD entry: KCONAZ) (Peeters et al., 1979) and rifaximin (CSD entry: ZEFMUU) (Braga et al., 2012) and then searched for these in our ASDs.To quantify similarities in the molecular placements, we introduced vectors in each API defined by two atoms that we deemed unlikely to change conformation within the molecule or in interactions with surrounding ones.As cyclohexane rings have a rigid form, and therefore unlikely to bend, we used two opposite carbon atoms in these rings in every drug (see Fig. 1).In crystalline structures, molecules are strictly organized, with a limited number of distances between the same segments in the neighboring molecules and angles between them.This is not the case in ASDs.By testing the distances between the vectors and their collinearity, we could use MD to evaluate the similarity of API clusters to those of the crystalline drugs.
To characterize the crystalline structures, supercells of 3 × 3 × 3 repetitive units were made in PyMOL software (Delano, 2002) using the .ciffiles downloaded from the CCDC and sets of commands shared by Thomas Holder.(Holder, 2010) The x-, y-and z-components were evaluated for distances between the centers of geometry of the vector-containing ring (Fig. 1) in the central molecule and those in the surrounding molecules.In addition, the angles between the in-ring vectors, as well as between the vectors normal to the selected benzene rings, were analyzed with a python script (depicted with purple and orange arrows in Fig. 1).These characteristics, translational (Δx, Δy, Δz) and rotational (two angles between the vector pairs, were used as references for analyzing the MD trajectories in a MATLAB script.To be considered a pre-crystalline nucleus, the molecule had to be found at the correct place with respect to the reference molecule, assume a characteristic angle towards another, and be rotated correctly (as tested with the normal vectors).If a third molecule also met the criteria in relation to any molecule already part of a cluster, it was added to that cluster.
We introduced three parameters to control the tolerance in the clustering process.One was the maximal search distance within which the clustering would be allowed (it was set to be 1 nm).The second was an allowed deviation in the characteristic displacements (here set to 0.2 nm.The third parameter was the angle discrepancy (here chosen to be 5 • for both in-ring and normal vectors).With help of the first parameter, we could control the search range.For instance, if a perfect alignment was found on a relatively long distance, but the organization of molecules in between the aligned ones did not match a crystalline structure, it would not make sense to look for alignment on such a long scale.The second controlled precision of the relative placement of the two molecules.Finally, the angle control filtered out "incorrectly" rotated molecules.Periodic boundary conditions were taken into account in the MATLAB script.

Amorphous solid dispersion characterization
Fig. 2 and Tables S2-S6 give an overview of the solid-state form, T g event, and value of freshly prepared solid dispersions (the latter in supplementary tables only) with different drug loadings, as observed with DSC.Rifaximin could be transformed into an amorphous form with the aid of BLG regardless of the blend ratio (Fig. S1a, Table S2).All samples showed only one T g event, and the T g value of rifaximin-BLG ASD decreased with increasing drug loading.
In the furosemide-BLG solid dispersions, all formulations showed the characteristic amorphous halo with drug loadings from 10 % to 90 % (see XRPD diffractograms in Figure S1b).Formulations showed a single T g for the drug loadings in range from 10 % to 40 %, indicating that they were homogeneous amorphous single-phase systems.In contrast, the formulations with a drug loading over 40 % showed a second T g event, suggesting that these were heterogeneous amorphous mixtures (Fig. 2b and Table S3).At drug loadings above the 40-50 % threshold, a measurable drug-rich cluster formed, corresponding to a second T g .The lower T g could be attributed to the drug-rich amorphous phase, and the higher T g to a drug-saturated mixture of drug-BLG.
In the ibuprofen-BLG solid dispersions, samples with drug loadings from 10 % to 50 % showed the characteristic amorphous halo (Fig. S1c), suggesting a complete amorphization of the drug.However, for those exceeding 50 %, crystalline peaks of ibuprofen were detected (Fig. S1c), indicating that the complete transformation of ibuprofen into an amorphous form was not achieved with the aid of BLG.Samples with drug loading below 30 % showed a single T g event, indicating that they were homogeneous amorphous single-phase systems (Fig. 2c and Table S4).In contrast, two T g events were observed for ibuprofen-BLG ASD at 40 and 50 % drug loadings, indicating that the samples were heterogeneous amorphous mixtures.Moreover, ibuprofen-BLG ASD at 50 % drug loading exhibited a melting point in the DSC curve, suggesting the occurrence of recrystallization before reaching the melting point temperature.Additionally, crystallization peaks were also detected by XRPD at 50 % drug loading after one day at room temperature, further indicating that the sample at that loading was not stable and had a high propensity for recrystallization.Samples with drug loading exceeding 50 % could not achieve complete amorphization, likely due to the recrystallization of the drug-rich cluster phase during the spraydrying process.
Griseofulvin could be transformed into an amorphous form with the aid of BLG at loadings from 10 % to 80 % (Figure S1d).Samples with drug loading below 40 % showed only one T g event, and the values decreased with drug loading, indicating that they were homogeneous amorphous single-phase systems.Conversely, two T g events were detected for samples with drug loading between 40% and 80%, suggesting that they were heterogeneous amorphous mixtures (see in Fig. 2d or Table S5).The lower T g could be attributed to the drug-rich amorphous phase, while the higher T g could be assigned to a drugsaturated mixture of drug-BLG.Interestingly, the lower T g values of griseofulvin-BLG ASDs closely resembled the T g of pure amorphous griseofulvin (shown with a dashed line in Fig. 2d), whereas the lower T g of furosemide-BLG ASD was clearly higher than the T g of pure amorphous furosemide.
Like griseofulvin, ketoconazole also features only hydrogen bond acceptors.In ketoconazole-BLG solid dispersions, all formulations showed the characteristic amorphous halo with drug loading from 10 % to 90 % (Fig. S1e).Samples with drug loadings from 40 % to 90 % showed two T g events, suggesting that they were heterogeneous amorphous mixtures (Fig. 2e and Table S6).The lower T g could be attributed to the drug-rich amorphous phase, while the higher T g could be assigned to a drug-saturated mixture of drug-BLG.For samples with 10 % to 30 % drug loadings, only one T g event was observed.However, the T g value of ketoconazole-BLG ASD at 30 % drug loading exhibited a significant decrease compared to samples with 10 and 20 % drug loadings, approaching the T g value of pure amorphous ketoconazole.

Physical stability
Rifaximin-BLG ASD with all drug loadings remained stable at the end of the experiment (Table 2).Previous research indicates that pure amorphous rifaximin (without stabilizing agent) exhibits remarkable physical stability-it can maintain stability for at least three months under 75 % relative humidity conditions at 40 • C (Viscomi et al., 2016).
Furosemide-BLG ASDs remained stable at drug loadings of up to % for at least 12 months under dry conditions at 40 • C. The recrystallization was detected for furosemide-BLG ASDs with drug loading of % and higher after 6 months.These findings are also consistent with the physical stability results of indomethacin-BLG ASD (Kabedev et al., 2022).Ibuprofen-BLG ASDs remained stable at drug loadings up to 40 % for at least 12 months.However, recrystallization was detected for ibuprofen-BLG ASD at 50 % drug loading during the physical stability experiment.At drug loadings higher than 50 %, crystalline phase was present already in the freshly prepared samples.Recrystallization of griseofulvin in the ASDs was detected for samples with drug loadings of 60 % and higher after 6 months.Ketoconazole-BLG ASDs with different drug loadings remained stable for at least 12 months at 40 • C under dry conditions.In contrast, amorphous ketoconazole, without any carrier, recrystallized within 3 days following storage at 40 • C. Remarkably, even at low BLG concentrations, BLG exhibited a significant ability to inhibit ketoconazole recrystallization.

Molecular dynamic simulations
In this section, we present the results from several analyses performed on all five drug molecules.

Diffusivity of the drug molecules
One of the most commonly used metrics to analyze stability of molecules with MD simulations is molecular drug diffusivity.High average diffusivity of the API, 〈D〉, might indicate the ability of the molecules in ASD to rearrange on longer timescales than simulations, which can result in re-crystallization.
The diffusivity averaged over all API molecules was the highest at the maximal drug loading (Table 3).Rifaximin molecules had the lowest values among all the molecules.Even at 90 % drug loading, its average diffusivity did not exceed 1 × 10 − 10 cm 2 /s.The highest 〈D〉 values for each drug loading were in the systems with ibuprofen and ketoconazole,  X. Zhuo et al. with ibuprofen always being the more mobile of the two; it had 27 times higher mobility than rifaximin.A comparison of the diffusivity values and patterns with the experimentally observed stability of the drugs showed that a value from MD simulations of approximately 1 × 10 − 10 cm 2 /s is a reasonable threshold of formulation instability.This particular cutoff results in the least number of mismatches between computational and experimental analysis for the set of APIs in this study.
Nevertheless, the molecular diffusivity for all drug loadings of ketoconazole, as well as griseofulvin at 60 % and ibuprofen at 40 %, does not match this criterion when compared to the experimental stability.The ketoconazole values are all greater than the cutoff, whereas experimentally the ketoconazole-BLG ASD remained stable for 12 months even at the highest drug loading.Recrystallization of griseofulvin at 60 % drug loading was observed already after six months; thus, the relatively low diffusivity value of 0.67 × 10 − 10 cm 2 /s does not match the experiment.Finally, at 40 % drug loading, ibuprofen was stable in the experiment, whereas in the simulations 〈D〉 was higher than 1 × 10 − 10 cm 2 /s.
It is important to note that because of the wide range of diffusivity values of different molecules in the system (Fig. 3a), the standard deviations of the diffusivity values are large.Especially in the ibuprofen-BLG and ketoconazole-BLG systems, the portion of the highly mobile molecules was substantial.For furosemide, griseofulvin and rifaximin ASDs, more than three-quarters of the molecules were very slow, as the leftmost points in panel (a) suggest.At the same time, the border between the closest 50 % and the most remote 50 % of the APIs is beyond one nm for all drugs (Fig. 3b).This clearly indicates-at least for these three APIs-that relatively slow drug molecules are not only in the first drug shell around the proteins (API molecules with <1 nm distance from the BLG protein surface).Moreover, on average there was very little correlation between the diffusivity of the individual molecules and their respective distances from the BLG surface.These findings served as a motivation to look in more detail at the distribution of the relatively more mobile drug molecules across the 90 % drug loading simulations.
Starting with visual analysis (Fig. 4), we differentiated the molecules to make APIs with low diffusivity more transparent.This allowed us to investigate in detail whether the most mobile molecules were predominantly distant from the BLG and whether any clear patterns could be observed.As the visual analysis showed, the distribution of opaque molecules was rather even across the ASDs, covering both the first drug shell (depicted in blue in Fig. 4) and more remote regions of the simulation boxes (red in Fig. 4).Moreover, the mixed nature of the low and high diffusivity API molecules, at different distances, was characteristic for all of the drug models.Together with a near-zero correlation coefficient between diffusivity and distance from the protein, it indicates the complexity of the stabilization pattern.

Hydrogen bonds
The next step to the investigation of the diversity among diffusivity values was to analyze the hydrogen bonds within the ASDs.The presence of hydrogen bond networks was previously found to be a crucial factor for indomethacin-BLG ASDs (Kabedev et al., 2022).
Furosemide had an average of 1.78 and 2.05 hydrogen bonds with proteins at drug loadings of 40 and 60 %, respectively (Table 4).It has both donors and acceptors and the molecule is relatively small, which a R: recrystallization.b C: crystallinity found in freshly prepared samples.c A: amorphous.
Abbreviations: FURfurosemide, GRIgriseofulvin, IBUibuprofen, KETketoconazole, RIFrifaximin, DLdrug loading.increases the chance to form a hydrogen bond with a BLG protein.At the drug loading of 90 %, the number dropped to 0.16, as most of the molecules are not in direct contact with the stabilizing protein.Rifaximin followed a similar trend.At drug loadings of 40 and 60 % it had 0.97 and 1.13 hydrogen bonds with the protein, and 0.24 at 90 %.Rifaximin also has a number of donors and acceptors, but the molecule is rather big (Table 1).Among the studied APIs, griseofulvin and ketoconazole had the lowest average number of hydrogen bonds with BLG.These are molecules with only hydrogen bond acceptors (six in each case).Ibuprofen had an intermediate number of hydrogen bonds with the proteins, starting with 0.89 at the lowest drug loading and 0.08 at the highest drug loading.Thus, most of the ibuprofen molecules might be stabilized, among other factors, by hydrogen bonds at the 40 % drug loading, whereas a relatively minor part of the ibuprofen molecules would be in contact with the surface at 90 %.
Another potentially important factor for ASD stabilization is the formation of hydrogen bonds between the drug molecules themselves (Table 5).Clearly, there are no such bonds between the APIs only carrying acceptors (griseofulvin and ketoconazole).Furosemide had on average more than one API-API hydrogen bond at all drug loadings, reaching a maximum of 2.98 at 40 % drug loading.Despite the fact that rifaximin has 5 donors and 12 acceptors, the average number of hydrogen bonds formed between these APIs was lower than that for furosemide.Finally, ibuprofen had an average number of hydrogen bonds per molecule below 1 for all drug loadings.
We found no significant correlation between the number of hydrogen bonds formed by individual molecules and the diffusivity of that molecule.The correlation coefficients were close to zero for all drugs, both for API-BLG and API-API hydrogen bonds.

API interactions with BLGs
The relative affinity to the proteins among five APIs was estimated via the determination of the number of protein-protein, protein-API, and protein-water contacts.In case contacts with APIs were underrepresented as compared to other drugs, we assumed lower protein-drug affinity.The number of BLG-BLG contacts is presented as the percentage of the total number of contacts between the BLGs and all surrounding molecules (Fig. 5).
Griseofulvin and ketoconazole were the drugs that mixed least with the protein, with approximately 1.26 and 1.22 times fewer number of fractions of BLG-API contacts than the other drugs, respectively.At the same time, the percentage of the protein-protein interactions for these two APIs was the highest at all drug loadings.For instance, at the drug loading of 90 %, inter-BLG contacts were 10.61 and 12.54 % for griseofulvin and ketoconazole, whereas for furosemide and rifaximin the values were 2.92 and 2.74 %.Even though ibuprofen had an intermediate value of 6.36 % of all contacts being BLG-BLG at 90 %, the drug was the most miscible with BLG, as can be seen at the longest pale orange bars in Fig. 5.The differences obtained here imply differences in the affinity of the drugs to the proteins, something that in turn can affect the ability of the proteins to stabilize drug molecules.These data are in good agreement with the number of hydrogen bonds between ketoconazole/griseofulvin and BLGs, as shown in Table 4.When low affinity is estimated this way, the number of bonds formed with the protein is also relatively low.

Crystalline structure patterns observed in silico
As diffusivity and hydrogen bond networks might not fully explain the stability of the drugs in BLG ASDs, we introduced a crystalline pattern search through the simulated ASD systems.
The crystalline structures from the CCDC were used as references for the vectors between neighboring molecules.In cases where the distances between the molecules in all three dimensions and the angles (as defined in Fig. 1) matched those in the crystal within the tolerance ranges, the molecules were considered part of a possible pre-crystallization cluster.Fig. 6 (left panel) shows such molecules from a 90 % drug loading in the griseofulvin-BLG ASD system, together with the proteins.The right panel is a more schematic representation of the vectors depicting the alignment of the molecules in the clusters.
A comparison of the identified API clusters for all five drugs showed Fig. 4. Distribution of the mobile molecules across the 90 % drug loading simulation boxes.The molecules with 〈D〉 of 5 × 10 − 10 cm 2 /s and lower are completely transparent, whereas those with diffusivity of 5 × 10 − 9 cm 2 /s and higher are fully opaque.The range is chosen in such a way that the major part with the slowest molecules present in all systems is omitted, while there is a common upper border for all systems (more detail in Section 2.3.3.1.).Blue indicates the molecules in contact with BLGs and red the ones without such contacts.These are relatively fast molecules both in the vicinity, and at a greater distance, from the BLGs in all four panels.This applies to ketoconazole and ibuprofen as well, but due to the opaqueness of many mobile molecules, the distribution of the transparent molecules is obscured.The box for ibuprofen is not shown, as it is fully opaque red and carries no information.Abbreviations are the same as in Fig. 3. Abbreviations are the same as in Table 2. Abbreviations are the same as in Table 2.
X. Zhuo et al. that a relatively low number of pre-crystalline clusters were present in the furosemide-BLG ASD at 90 % drug loading (Fig. 7).This is in line with the recrystallization observed experimentally at the high drug loadings.Also, in agreement with XRPD data, griseofulvin showed an even higher number of pre-crystalline clusters.Ibuprofen exhibited the highest number of clusters, matching with the lowest stability observed among all drugs.Ketoconazole and rifaximin are the APIs that remained amorphous in the experimental studies, and would presumably have fewer clusters.Ketoconazole indeed had the fewest number of clusters among the drugs, but eventually clusters were formed throughout the MD trajectory for both ketoconazole-BLG and rifaximin-BLG ASDs (Fig. 8).In summary, the number of pre-crystalline clusters in the simulation boxes was in good agreement for all five drugs with the experimental studies on the physical stability of the ASDs at 90 % drug loading.
Another potential in silico criterion for the evaluation of experimental stability is the average number of molecules present in different clusters throughout the simulations.Here, only griseofulvin-BLG and ibuprofen-BLG ASDs contained clusters with at least n = 3 molecules at 90 and 60 % drug loadings (Fig. 8).The ibuprofen systems had a significantly greater number of clusters in all cases.No cluster of more than two molecules was found in furosemide-BLG system, but the average number of molecules in clusters of 2 was an order of magnitude greater than in ketoconazole-BLG and rifaximin-BLG ASDs.

Discussion
The experimental results for furosemide were consistent with the outcome of our previous study on indomethacin-BLG ASD stability.Specifically, at drug loadings up to 40-50 %, indomethacin and BLG Fig. 5. Percentage of number of contacts made by BLGs with other protein molecules, APIs, water, and ions in different API-ASD systems.The BLG proteins had the largest fraction of contacts with molecules like itself in the griseofulvin-BLG and ketoconazole-BLG ASDs, while the respective fractions of contacts with APIs were the lowest.This is an indicator of the low affinity of these two drugs to the BLG.Abbreviations are the same as in Fig. 3. X. Zhuo et al. formed a homogeneous amorphous single-phase system with a single T g .The ASD also remained stable for 12 months even at a drug loading of 70 %.It is tempting to conclude that stabilization is explained by low molecular mobility in the first drug shell-at least for systems with drug loadings less than 50%, like in our previous study on indomethacin.
However, the situation is more complex than this, as we discuss in the following for each drug.
The stability mechanism of ibuprofen-BLG ASDs with a drug loading below 30% (homogeneous amorphous system) could be primarily ascribed to the limited mobility of drug molecules within the first drug Fig. 7. Possible crystallization sites (clusters of molecules) across the ASDs with five drugs at 90 % drug loading, performed at the final frames of the respective MD trajectories.The analysis was done with a cut-off distance of 1.0 nm (see methods) and with angle and displacement thresholds of 5 • and 0.2 nm, respectively.The clusters are distinguishable by unique combinations of color-coding and marker symbols for each group.Abbreviations are the same as in Fig. 3. Fig. 8. Averaged cluster size distribution for three simulated drug loadings.Cluster size of two with an average number of molecules equal to one indicates that one cluster was present in half of the MD trajectory (as, for example, for furosemide at 90 % DL).Analysis was done over the last 150 ns of the simulations.Abbreviations are the same as in Fig. 3. shell.But for the sample with 40 % drug loading, a second T g appears.The stabilization mechanism here could be a combination of steric confinement provided by the BLG matrix and hydrogen bonding.Recrystallization was detected for ibuprofen-BLG ASD at 50% drug loading during the physical stability experiment.This higher propensity of the ibuprofen in drug-rich clusters to undergo recrystallization is primarily due to the extremely low T g of pure ibuprofen.This leads to the relatively high mobility of ibuprofen molecules, also confirmed with computational data (Table 3).Given that ibuprofen only has one functional group that can form hydrogen bonds with BLG (Figure S4), it probably cannot form a hydrogen bond network that extends to molecules in the second drug shell and beyond.Even though the carboxylic acid group in ibuprofen in theory can act as a hydrogen bond donor, it might be sterically hindered to act as an acceptor once it is involved in hydrogen bonding with BLG.Hence, the molecules not present in the first drug shell are not likely to be stabilized by a hydrogen bond network.This leaves them prone to recrystallization, due to the low Tg of ibuprofen and the poor physical stability of pure amorphous ibuprofen.
Even though rifaximin has multiple donors and acceptors, their placement around the molecule does not seem to facilitate the formation of a far-extending hydrogen bond network.The computational analysis showed that the average number of hydrogen bonds formed by rifaximin molecules with each other is lower than that of furosemide and does not even reach 1 at 90% drug loading.At the same time, we know that pure amorphous rifaximin can also stay stable for at least three months.It brings us to the conclusion that a great deal of external energy is required to rearrange the molecular organization from an amorphous one to a crystalline one.Factors associated with this high energy cost include: Van der Waals forces and miscibility of the API with the stabilizing agent (Walden et al., 2021), ionic interactions (Maniruzzaman et al., 2013), rotatable bonds within the drug molecules (Edueng et al., 2021).In addition, another key feature of rifaximin, which likely determines its high stability at all drug loadings is the high T g and T m .
The difference in the T g values of griseofulvin-BLG and ketoconazole-BLG ASDs is small compared to the pure amorphous griseofulvin and ketoconazole (Fig. 2d and e).This contrasts with the significantly higher T g value of furosemide-BLG ASD compared to pure furosemide (Fig. 2b).The difference might originate from the presence of hydrogen bonds between furosemide molecules.Since griseofulvin and ketoconazole only possess hydrogen bond acceptors, the hydrogen bond network originating from the BLG cannot form in the drug-cluster phase of these ASDs, resulting in T g values similar to the pure amorphous griseofulvin.The poor stability of the griseofulvin-BLG samples at drug loadings of 60 % and higher can therefore be explained, at least partially, by the lack of a hydrogen bond network.The samples at drug loadings of 40 and 50 %, were not homogeneous, but the amorphous state of griseofulvin remained stable for at least 12 months.This may have been a result of the steric confinement provided by the BLG matrix.Interestingly, griseofulvin and ketoconazole, which only have acceptors, could potentially at least have been stabilized in the first drug shell, as BLG has numerous hydrogen bond donors.However, their respective average numbers of hydrogen bonds with the proteins are 0.46 and 0.62 at the minimal drug loading (Table 4).This can, to some extent, be explained by the low affinity of the drugs for the protein (Fig. 5).
Interestingly, the lower T g of ketoconazole-BLG samples with 50 % to 90 % drug loadings was slightly lower than that of pure amorphous ketoconazole (~44 • C).Similar results were found in previous studies where the addition of PVP (10 % w/w) or Eudragit E PO (25 % w/w) to ketoconazole also lowered the T g compared to the T g of pure amorphous ketoconazole (Ullrich and Schiffter, 2018).The stability of the ketoconazole-BLG ASD system at high drug loadings is unexpected considering that ketoconazole only features hydrogen bond acceptors.This implies that the hydrogen bond network for stabilizing the amorphous ketoconazole does not form in the drug-rich regions.The underlying stabilization mechanism might therefore be attributed to other types of interactions.Previous studies showed that polyacrylic acid (PAA) is also an effective polymer for inhibiting recrystallization of ketoconazole, even at low polymer content.Ketoconazole-PAA ASDs with 4 % PAA remained amorphous for >1 year after storage at 40 • C under dry conditions (Mistry et al., 2017), and ketoconazole-PAA ASDs with 30 % PAA remained amorphous for > 6 months after storage at 40 • C/75 % RH. (Sarpal et al., 2020) PAA seems to inhibit the recrystallization of ketoconazole by reducing molecular mobility, mainly via ionic interactions between the polymer and drug.In the ketoconazole-BLG ASDs, the presence of diverse amino acids within BLG might also lead to the formation of structural units that interact with ketoconazole through non-hydrogen bond electrostatic interactions.(Fung et al., 2018) Apart from that, the formation of hydrogen bonds between griseofulvin/ketoconazole and residues on BLG or water can induce changes in electron cloud density in the structure of the drugs.This potentially leads to the formation of weak, induced hydrogen bonds in the molecules (Aitipamula et al., 2012;Chen et al., 2020) and thereby stabilizes these ASDs.Less likely, but worth considering, is that residual solvent such as acetic acid from the spray drying remains in the ketoconazole-BLG ASD.This might also facilitate the formation of a stable co-amorphous phase.Finally, all the molecules in this study were simulated in their un-ionized form, which might not be an accurate representation of an acidic ketoconazole.
Analysis of the hydrogen bonds reported in Section 3.3.2demonstrated that the average number of hydrogen bonds was rather low at high drug loadings for four of the five molecules.Only furosemide had an average number of API-API hydrogen bonds higher than one at 90% drug loading (Table 5).Similar to indomethacin, furosemide has a rather broad distribution of donors and acceptors across a relatively small molecule (Fig. S4).The average number of hydrogen bonds with BLGs was also the highest for furosemide, whereas for the other four molecules it was lower than one even at a drug loading of 40 %.Thus, for samples with drug loadings ranging from 50 % to 70 %, the stabilization mechanism could be the combination of the hydrogen bond network and the steric confinement of the molecules.
The average diffusivity of the molecules is another criterion used to computationally evaluate the stability of ASDs.In Section 3.3.1 we analyzed this and found a certain match in patterns of the values with the physical stability data.We chose the threshold of stability to be 1 × 10 − 10 cm 2 /s to get the best fit with the experimental data.Rifaximin remained stable after 12 months at all drug loadings, and in simulations, the mobility of the molecules was below the reference value at all drug loadings (Table 3).Furosemide also showed a perfect match with experimental data, as it becomes unstable at loadings higher than 70 %.Ketoconazole was much more mobile in the simulations than expected-in the stability experiments, it remained amorphous at all loadings after one year.Finally, ibuprofen and griseofulvin showed mediocre agreement with the experimental data.The stabilities of ibuprofen at 40 % and griseofulvin at 60 % drug loading, were underestimated, at least when only looking at average diffusivity.Clearly, diffusivity is a useful tool for a quick assessment of ASD stability, but alone it might be not sufficient.There is no clear general rule behind the selection of a threshold value for a robust stability prediction of new compounds except for determining molecular diffusion for numerous drug models and analyzing the data (1 × 10 − 10 cm 2 /s worked best for this study).Even our best fit did not lead to a full match with the stability data.In other words, high mobility in the simulations is not necessarily a proxy for low stability of the molecules in an experimental setup.
MD simulations have multiple limitations.For instance, most force fields cannot capture covalent bond formation or subtle induced hydrogen bonds.As discussed above, ketoconazole and griseofulvin could theoretically form induced hydrogen bonds, which might be an important factor for the stabilization of those ASDs.Higher diffusivity of individual molecules on the nano-to micro time and length scale might not necessarily be extrapolated to the instability of a one-year storage process in the laboratory.For example, drug molecules mobile on the X.Zhuo et al. small scale can be "trapped" by the systematic arrangement of BLGs.As a result, their diffusivity would not represent low stability.Finally, there might not be a trivial connection between higher diffusivity and lower stability of the simulated ASD.On the one hand, it is intuitive that more mobile molecules can align themselves faster and form precrystallization nuclei.Some mechanical models also confirm this relation (Ojo and Lee, 2021).On the other hand, low mobility might indicate that the API molecules are probably already in a good position to form crystalline structures and that these molecules crystallize.A more complex analysis protocol, including crystallinity pattern check, is needed to reach more confident conclusions.
In the last in silico part of this study, we analyzed molecules in the ASDs for the presence of pre-crystalline nuclei.Characteristic distances between the drugs, as well as the angles between certain intramolecular vectors, were extracted from the Cambridge Crystallographic Data center database.We then searched for such patterns in the ASD simulations boxes and found the results to be in good agreement with the physical stability experiments.Ketoconazole and rifaximin had the lowest number of pre-crystalline clusters with aligned molecules, indicating that there were few crystallization nuclei, if any (Figs.7 and 8).Furosemide had a higher number of aligned molecules at the characteristic distances from each other (Fig. 7).Griseofulvin molecules formed significantly more pre-crystalline nuclei across the box.It also had clusters of more than two molecules at 60% and 90% drug loadings (Fig. 8).This can be interpreted as a higher probability of crystal formation from the current structure.Ibuprofen showed the highest number of pre-crystalline clusters across the box.It might be closely related to its high diffusivity, as the molecules are the most mobile and can more easily rearrange into energetically more favorable structures.This also puts the ibuprofen simulation in good agreement with the experimental data.The ibuprofen-BLG ASD was the least stable, not even obtainable in an amorphous form at 90% drug loading, and it was not stable for 12 months at 60% drug loading.

Conclusion
In this study, we have examined ASDs of five drugs with BLG at multiple drug loadings.XRPD data showed that ketoconazole and rifaximin remained amorphous at least for 12 months at drug loadings of up to 90 % (of total weight).Furosemide remained stable at loadings of up to 70 %, griseofulvin up to 50 %, and ibuprofen only at 40 % drug loading.We modeled the ASDs with molecular dynamics and evaluated the diffusivity of the drug molecules, API-API and API-BLG hydrogen bonds, and affinity of APIs to BLG.Finally, we introduced a search strategy for crystalline patterns through the simulation box to find precrystalline nuclei.Using a combination of experimental techniques and computer simulations, we showed that hydrogen bond networks are a stabilizing mechanism in BLG-based ASDs.This was true for drugs with either both donors and acceptors as with furosemide or indomethacin, and for drugs with acceptors and induced donors, as can potentially be the case for ketoconazole or griseofulvin.However, relative placement of the donor and acceptor sites should be broad enough across the molecule to allow formation of the network.Other factors affecting the ability to remain amorphous are steric confinement, high T g , and the geometry of the drug.
The combined computational approach, based on (1) analysis of API diffusivity, (2) analysis of the hydrogen bonds (specifically hydrogen bond networks), and (3) crystallinity pattern searches, has proven to be a powerful tool.It potentially can be further developed for in silico ASD stability predictions, by taking into account the T g and T m of the API as well as the drug's ability to form induced hydrogen bonds and ionic bonds.

Fig. 1 .
Fig. 1.Drug molecules with vectors used to define alignment and match with patterns from the crystalline structures.Purple arrows: in-ring vectors.Orange arrows: vectors normal to the rings.Abbreviations: FUR -furosemide, GRIgriseofulvin, IBUibuprofen, KETketoconazole, RIFrifaximin.

Fig. 2 .
Fig. 2. Experimentally obtained T g values of rifaximin (RIF), furosemide (FUR), ibuprofen (IBU), griseofulvin (GRI) and ketoconazole (KET) ASDs as a function of drug loading (DL).Dashed lines show the pure drug T g values.At the drug loadings with two T g values, squares mark the higher of them.

Fig. 3 .
Fig. 3. Distributions of individual API molecule diffusivities D (a), and distance from the BLG protein surface d (b), measured at 90 % drug loading.Dotted lines in panel (b) depict the average distance splitting half of the API population closer to the BLGs from those more remote from it.Abbreviations: FURfurosemide, GRIgriseofulvin, IBUibuprofen, KETketoconazole, RIFrifaximin.

Fig. 6 .
Fig. 6.Visualization of the crystallinity patterns in the griseofulvin-BLG ASD at 90 % drug loading.In the left panel, the snapshot from the simulation shows only molecules meeting the relative placement and cosine similarity criteria (see methods).The BLG proteins are shown in light orange.The right panel shows only the vectors (see Fig. 1) from the drug molecules in the left panel.

Table 1
Physicochemical properties of the model drugs.

Table 2
Physical Stability of ASDs samples at 40 • C, dry conditions Data are valid for both 6 and 12 months in all cells.

Table 3
Average diffusivity (<D>, × 10 − 10 cm 2 /s) of API molecules in ASDs as observed from MD simulations.Values that do not agree well with the observed physical stability trends are shown in boldface.Abbreviations are the same as in Table2.

Table 4
Average number of hydrogen bonds formed by an API with BLG.

Table 5
Average number of hydrogen bonds formed by an API with other APIs.