pharmaceutics Article Rational Design of Hydrogels for Cationic Antimicrobial Peptide Delivery: A Molecular Modeling Approach Alfredo Pereira 1, *, Elizabeth Valdés-Muñoz 2 , Adolfo Marican 3,4,5 , Gustavo Cabrera-Barjas 6 , Sekar Vijayakumar 7 , Oscar Valdés 8 , Diana Rafael 9,10,11 , Fernanda Andrade 9,10,12 , Paulina Abaca 4 , Daniel Bustos 8,13, * and Esteban F. Durán-Lara 4,5, * 1 Departamento de Química Orgánica y Fisicoquímica, Facultad de Ciencias Químicas y Farmacéuticas, Universidad de Chile, Santiago 8380544, Chile 2 Doctorado en Biotecnología Traslacional, Facultad de Ciencias Agrarias y Forestales, Escuela de Ingeniería en Biotecnología, Universidad Católica del Maule, Talca 3480094, Chile 3 Instituto de Química de Recursos Naturales, Universidad de Talca, Talca 3460000, Chile 4 Bio & Nano Materials Lab, Drug Delivery and Controlled Release, Departamento de Microbiología, Facultad de Ciencias de la Salud, Universidad de Talca, Talca 3460000, Chile 5 Center for Nanomedicine, Diagnostic & Drug Development (ND3), Universidad de Talca, Talca 3460000, Chile 6 Unidad de Desarrollo Tecnológico (UDT), Universidad de Concepción, Parque Industrial Coronel, Coronel 3349001, Chile 7 Marine College, Shandong University, Weihai 264209, China 8 Centro de Investigación de Estudios Avanzados del Maule (CIEAM), Vicerrectoría de Investigación y Postgrado Universidad Católica del Maule, Talca 3460000, Chile 9 Drug Delivery & Targeting, Vall d’Hebron Institut de Recerca (VHIR), Universitat Autònoma de Barcelona (UAB), 08035 Barcelona, Spain 10 Centro de Investigación Biomédica en Red de Bioingeniería, Biomateriales y Nanomedicina (CIBER-BBN), Instituto de Salud Carlos III, 28029 Madrid, Spain 11 Functional Validation & Preclinical Research (FVPR)/U20 ICTS Nanbiosis, Vall d’Hebron Institut de Recerca (VHIR), Universitat Autònoma de Barcelona (UAB), 08035 Barcelona, Spain 12 Departament de Farmàcia i Tecnologia Farmacèutica i Fisicoquímica, Facultat de Farmàcia i Ciències de l’Alimentació, Universitat de Barcelona (UB), 08028 Barcelona, Spain 13 Laboratorio de Bioinformática y Química Computacional, Departamento de Medicina Traslacional, Citation: Pereira, A.; Valdés-Muñoz, E.; Marican, A.; Cabrera-Barjas, G.; Facultad de Medicina, Universidad Católica del Maule, Talca 3480094, Chile Vijayakumar, S.; Valdés, O.; Rafael, D.; * Correspondence: apereira.bioinfo@gmail.com (A.P.); dbustos@ucm.cl (D.B.); eduran@utalca.cl (E.F.D.-L.) Andrade, F.; Abaca, P.; Bustos, D.; et al. Rational Design of Hydrogels for Abstract: In light of the growing bacterial resistance to antibiotics and in the absence of the develop- Cationic Antimicrobial Peptide ment of new antimicrobial agents, numerous antimicrobial delivery systems over the past decades Delivery: A Molecular Modeling have been developed with the aim to provide new alternatives to the antimicrobial treatment of Approach. Pharmaceutics 2023, 15, infections. However, there are few studies that focus on the development of a rational design that 474. https://doi.org/10.3390/ is accurate based on a set of theoretical-computational methods that permit the prediction and the pharmaceutics15020474 understanding of hydrogels regarding their interaction with cationic antimicrobial peptides (cAMPs) Academic Editor: Chaoliang He as potential sustained and localized delivery nanoplatforms of cAMP. To this aim, we employed docking and Molecular Dynamics simulations (MDs) that allowed us to propose a rational selec- Received: 3 January 2023 tion of hydrogel candidates based on the propensity to form intermolecular interactions with two Revised: 19 January 2023 types of cAMPs (MP-L and NCP-3a). For the design of the hydrogels, specific building blocks were Accepted: 28 January 2023 Published: 31 January 2023 considered, named monomers (MN), co-monomers (CM), and cross-linkers (CL). These building blocks were ranked by considering the interaction with two peptides (MP-L and NCP-3a) as receptors. The better proposed hydrogel candidates were composed of MN3-CM7-CL1 and MN4-CM5-CL1 termed HG1 and HG2, respectively. The results obtained by MDs show that the biggest differences Copyright: © 2023 by the authors. between the hydrogels are in the CM, where HG2 has two carboxylic acids that allow the forming of Licensee MDPI, Basel, Switzerland. greater amounts of hydrogen bonds (HBs) and salt bridges (SBs) with both cAMPs. Therefore, using This article is an open access article theoretical-computational methods allowed for the obtaining of the best virtual hydrogel candidates distributed under the terms and according to affinity with the specific cAMP. In conclusion, this study showed that HG2 is the better conditions of the Creative Commons candidate for future in vitro or in vivo experiments due to its possible capacity as a depot system and Attribution (CC BY) license (https:// its potential sustained and localized delivery system of cAMP. creativecommons.org/licenses/by/ 4.0/). Pharmaceutics 2023, 15, 474. https://doi.org/10.3390/pharmaceutics15020474 https://www.mdpi.com/journal/pharmaceutics Pharmaceutics 2023, 15, 474 2 of 14 Keywords: cationic antimicrobial peptides; hydrogels; molecular docking; molecular dynamics simulations; virtual-screening; antibiotic resistance 1. Introduction About two million fatalities were attributed to antimicrobial resistance (AMR) in 2020 alone, a number that grew due to the COVID-19 pandemic’s indiscriminate antibiotic usage [1]. By 2050, it is predicted that over 10 million people will have died as a direct result of antibiotic resistance [2]. Nosocomial pathogenic micro-organisms with rising degrees of antibiotic resistance are referred to by the abbreviation ESKAPE for the Gram- positive species Enterococcus faecium and S. aureus and the Gram-negative K. pneumoniae, A. baumannii, P. aeruginosa and Enterobacter species. ESKAPE pathogens tend to inhabit hospital environments, and they are the most common cause of fatal infections in immuno- compromised patients [3,4]. Despite a vital need for novel antimicrobial drugs, just a few new viable antibiotics are currently in clinical development. In this context, AMPs appear as promising candidates in the fight against multi-drug resistance bacterial infections due to their broad range of activities and low toxicity [5]. AMPs are small proteins that have antibacterial broad-spectrum and immunomodulatory effects on infections produced by micro-organisms such as bacteria, viruses, fungi, or parasites. These peptides can be char- acterized according to their sequence-structure composition, e.g., net-charge, secondary and tertiary architecture, solubility, electrostatic components, and source of origin (nat- ural or synthetic) [6]. Unlike antibiotics, which disrupt the synthesis of peptidoglycan in bacterial cell walls or inhibit the DNA gyrase, AMPs operate by creating membrane pores and attaching to molecules that are precursors to membranes [7,8]. One of the most abundant AMPs is cAMP, which can disrupt the bacterial membrane, resulting in cell death via osmotic shock [6]. The remarkable spectrum of functions that cAMP may perform includes antibacterial action against both Gram-positive and Gram-negative bacteria [9]. Furthermore, cAMPs are auspicious candidates to replace conventional antibiotics; how- ever, they possess some drawbacks, such as their short half-life due to enzyme degradation, low bioavailability and solubility, and low activity in physiological conditions, among others [10,11]. An effective action to resolve the problem is to apply drug delivery nano- systems (DDNSs) that allow a sustained or prolonged release of AMPs. Thus, a rational design of effective cAMP-delivery nano-systems that provide protection and an effective release is essential for the successful implementation of cAMP treatments [12]. Some studies have proved the efficacy of combining hydrogels and cAMPs, enhancing the antimicrobial potency compared to the free-cAMP alternative and low cell toxicity [13,14]. Hydrogels have received growing attention in recent years as DDNSs due to their exclusive properties, such as high biocompatibility, tunable release rate, and versatility with regard to being loaded with different molecules. Consequently, to rationally design hydrogels, taking into consideration the interaction affinity with cAMPs is fundamental for achieving high-release performance. For example, unpredictable slow releases of drugs from hydrogels were attributed to some associative interactions between the drug and the matrix, impairing the release rate [15]. Stability and aggregation are other common difficulties regarding DDNS preparation, i.e., when an amphiphilic peptide is encapsulated into a hydrophobic system that includes an aqueous layer, it could leak from the systems [16]. The rational design of the on-demand formulations must be accurate in order to avoid the usual DDNS limitations. Under this new approach, injectable and stimuli-responsive hydrogels seem to be good candidates to solve these concerns. These platforms are highly attractive biomaterials in drug delivery and tissue engineering [17]. This class of hydrogels responds to small varia- tions in environmental conditions, such as light, ionic strength, magnetic fields, pH, and temperature. Among these, the most commonly studied hydrogels are responsive either to pH or temperature [18–20]. Until now, a wide variety of stimuli-responsive hydrogels has been reported; among them is the hydrogel of poly(N-isopropylacrylamide) (PNIPAM), Pharmaceutics 2023, 15, 474 3 of 14 the gold standard polymer in the field of thermoresponsive hydrogels (TRH) [21]. PNI- PAM has a lower critical solution temperature of 32−34 ◦ C [22,23]. Given its exclusive hydrophilic–hydrophobic transition behavior, various TRH PNIPAM-based hydrogels have been created for the controlled delivery of hydrophilic protein drugs. A special advantage regarding these hydrogels is that the loading of hydrophilic proteins can be achieved in aqueous media at a low temperature, which diminishes the chance of the denaturation of proteins in a hostile environment [18]. By using theoretical-computational methods such as docking and MDs, it is possible to understand and predict the interaction between hydrogels and different bioactive compounds such as AMP [24]. With these analyses, it could be possible to adjust the specific interaction level of hydrogel-cAMP, modifying the hydrogel structure with different building blocks. The above could guarantee the potential use of hydrogels as a depot system and a localized and sustained release platform of cAMPs. Therefore, in this study, we employed a set of docking and MDs that allowed us to propose a rational selection of hydrogel candidates based on the propensity to form intermolecular interactions with specific cAMPs as an innovative approach to the rational design of hydrogels applied to the sustained and localized release of cAMPs. 2. Materials and Methods 2.1. AMP Selection After a rigorous search, two cAMPs were selected—MP-L (ID: DRAMP18623), and NCP-3a (ID: DRAMP18636)—as cAMP models from the database http://dramp.cpu- bioinfor.org (accessed on 12 March 2022). The requirements of selection were (i) having a net positive charge, (ii) showing antimicrobial activity against the largest number of targets from the ESKAPE pathogens, and (iii) containing a number equal to or lower than 20 amino acids. The cAMP properties are shown in Table 1. The secondary and tertiary structures of each cAMP were predicted with PEP-FOLD 3.5 server [25]. Both cAMPs were energy minimized by using the conjugate gradient algorithm in a vacuum; then, each cAMP was embedded into a box of SPC (single-point charge) water molecules and chloride ions were added to neutralize the systems. Subsequently, the systems were equilibrated through 25 ns of all-atom MDs in an isobaric-isothermal NPT ensemble at 1 atm and 298.15 K with Maestro/Schrödinger suite [26] and OPLS3e as force fields [27]. The endpoint of cAMP in each simulation was collected for further analysis. Table 1. Peptides properties. MP-L or P1 NCP-3a or P2 Source Synthetic construct Synthetic construct Family Derived from the peptide Mastoparan Derived from CTX-1 Sequence INLKILARLAKKIL KLIFILSKTIPAGKNLFYKI Length 14 residues 20 Biological Activity Anti-Gram positive, Anti-Gram negative, Antifungal • S. aureus (ATCC 43300, ATCC 22953) • E. hirae (ATCC 10541) • S. aureus (ATCC 25923) • S. agalactiae (ATCC 13813) • S. pyogenes (ATCC 19615) • E. coli (ATCC 25922) • L. ivanovii (Li 4pVS2) • P. aeruginosa (ATCC 27853) • P. aeruginosa (ATCC 27853) • A. baumannii (ornithological, cloacal) Target Organism • K. pneumoniae (ATCC 13883) • K. pneumoniae (herpetological, cloacal) • A. baumannii (ATCC 19606) • B. cepacia (ATCC 17759) • C. albicans (ATCC 90028) • M. catarrhalis (ATCC 25238) • C. parapsilosis (ATCC 22019) • C. albicans (ATCC 10231) • C. glabrata (ATCC 90030) • M. pachydermatis (DSMZ 6172) Pharmaceutics 2023, 14, x FOR PEER REVIEW 4 of 14 Pharmaceutics 2023, 15, 474 4 of 14 ● C. albicans (ATCC 10231) Table 1. Cont. ● C. glabrata (ATCC 90030) MP-L or P1 ●NCP-3a M. or pachydermatis P2 (DSMZ 6172) Structure Structure predicted predicted The non-polar, polar and positively charged amino acids are depicted in gray, orange and blue, The non-polar, polar and positively charged amino acids are depicted in gray, orange and blue, respectively. The respectively. The bacteria highlighted in bold belong to the ESKAPE pathogen group. bacteria highlighted in bold belong to the ESKAPE pathogen group. 2.2. 2.2. Hydrogel Hydrogel Building Blocks Selection Building Blocks Selection The The selection of the selection of the building building blocks blocks such such as as MN, MN, CM, CM,andandCLCLwas wasmade madewith withthe theaim aim of evaluating the influence of the polarity and non-ionic interactions over the affinity of evaluating the influence of the polarity and non-ionic interactions over the affinity of of the hydrogel-peptide complex without negatively affecting the stimulus sensitivity the hydrogel-peptide complex without negatively affecting the stimulus sensitivity and and mechanical properties of mechanical properties of the the hydrogel. hydrogel. MN, MN, CM, CM, and andCLCLare areshown shownin inFigure Figure1.1. Figure Figure 1. 1. Chemical Chemical structures of MN, structures of MN, CM CM and and CL. CL. The structural units of the hydrogel systems, systems, called called “building “building blocks”, blocks”,were weresketched sketched in Maestro/Schrödinger Maestro/Schrödinger suite suite and and prepared prepared at pH 7 according according to to [24]. [24]. 2.3. Analysis 2.3. Analysis of of the the Affinity Affinity between between Hydrogel Hydrogel Building BuildingBlocks Blocksand andcAMPs cAMPs The building blocks were used as inputs jointly with eachcAMP The building blocks were used as inputs jointly with each cAMPto toperform performaavirtual virtual screening based on docking. For this, two different protocols and dockingprograms screening based on docking. For this, two different protocols and docking programswere were employed. First, employed. First,weweused usedthe theAutodock AutodockVina Vinaprogram program[28][28] toto dock dock each each MN, MN, CM, CM, andand CL CL against each cAMP (named Protocol 1). The Autodock Vina makes against each cAMP (named Protocol 1). The Autodock Vina makes use of a semi-flexible use of a semi- flexible approach approach in that itinonly thatallows it onlythe allows the(building ligand ligand (building blocks inblocks in this this case) to case) to beflexible. be made made flexible. In the second protocol (Protocol 2), the Induced-Fit docking (IFd) In the second protocol (Protocol 2), the Induced-Fit docking (IFd) approximation of glide approximation of glide software software was used was used with [29,30] [29,30] with the mosttherigorous most rigorous algorithm algorithm possible—XP possible—XP (Extra(Extra Preci- Precision)—and using the OPLS3e force field. In both protocols, the sion)—and using the OPLS3e force field. In both protocols, the final number of poses per final number of poses per pair (building block-cAMP) was set at 10. Then, all the poses through boxplots pair (building block-cAMP) was set at 10. Then, all the poses through boxplots comparing comparing the binding energy of each segment against the rest of its group in both protocols the binding energy of each segment against the rest of its group in both protocols were Pharmaceutics 2023, 15, 474 5 of 14 were statistically analyzed. The best building block belonging to each group (MN, CM, and CL) based on the binding affinity with cAMP were selected for the next step. 2.4. Ensemble and Simulations of Hydrogel-cAMP Complexes Four systems were ensembled: MN3-CM7-CL1 (termed as HG1) and MN4-CM5-CL1 (termed as HG2) hydrogels with each cAMP (MP-L and NCP-3a) independently. First, the LEAP program from AmberTools21 suite [31] was used to build substructures of hydrogels, considering the proportions of MN, CM and CL used in [24]. The HG1 substructure was constituted by two chains, each of 20 MN3 plus 5 CM7 (randomly distributed), cross-linked by a CL1. The HG2 substructure was constituted by 2 chains, each with 20 MN4 plus 5 CM5 (randomly distributed), cross-linked by a CL1. Then, the chains were equilibrated by MDs with the same parameters previously employed in the simulations of the cAMPs alone. The endpoint of each simulation was collected as input to build the hydrogel-cAMP complex. For this, Packmol software [32] was used to build a sphere of 60 Å of the radius where four hydrogel chains were placed, and another inner sphere of 50 Å of the radius where 10 cAMP molecules were positioned. Subsequently, the four systems were minimized, solvated, and neutralized with NaCl counter ions; finally, each system was simulated by 50 ns in the same condition as the previously described simulations. The MDs allow us to identify the main intermolecular interactions and forces between the hydrogels and cAMPs. All the analyses were performed over the entire simulation time and consisted of the gyration radius of each hydrogel, the secondary structure fluctuation of both peptides, the intermolecular contacts and non- bonding energies (van der Waals, electrostatic and total), and the interactions (HBs and SBs) between the hydrogels and the cAMPs. The intermolecular contacts were computed with an in-house TCL script and run in VMD software [33]. Regarding the non-bonding energies and interactions, these were calculated from the Maestro/Schrödinger suite and its Simulation Event Analysis tools. All the plots and statistical analysis were generated with R v.4.2.2 and RStudio v. 2022.07.0. 3. Results and Discussion 3.1. Ranking of Building Blocks by Molecular Docking Figures 2–4 present the building block rankings using the protocols based on the Autodock Vina (protocol 1) and IFd approximation of glide (protocol 2). The rankings of the MNs (Figure 2), CMs (Figure 3), and CLs (Figure 4) are shown considering the peptides MP-L and NCP-3a as receptors. For both strategies, the statistical comparison between building blocks were carried out by taking the 10 best poses of each docking run. Then, the selection of the best candidates was made based on the median value of the binding energy in both protocols individually. In both cases, the lower the value, the better the candidate. Overall, the results showed that the best candidates obtained coincide with the pep- tides studied; however, important differences were observed between the proposed strate- gies. For MN ranking, MN3 was the best candidate for both peptides using protocol 1, whereas MN4 was the best candidate for both peptides using protocol 2. For the CM ranking, using protocol 1, the best candidate for both peptides was CM7; however, in the case of protocol 2, the best candidate for both peptides was CM5. Finally, in the CL ranking, for both strategies and considering the two peptides studied, the best candidate turned out to be CL1. It is important to mention that when protocol 2 was used, CL4 also appeared as a good candidate, showing no significant difference with CL1; however, the median score value was still lower in CL1. The variations observed in the building blocks when comparing both protocols may be attributed to the specific docking method applied in each protocol. As mentioned above, AutoDock uses a semi-flexible variant; in contrast, IFd is based on the Induced-Fit hypothesis, which states that when a receptor (peptides) binds to a ligand (building blocks), it experiences conformational modifications, allowing for a better fit between the two molecules. The IFd method is thought to be more accurate than typical rigid-body docking methods because it takes into consideration the dynamic nature Pharmaceutics 2023, 14, x FOR PEER REVIEW 6 of 14 for a better fit between the two molecules. The IFd method is thought to be more accurate Pharmaceutics 2023, 15, 474 than typical rigid-body docking methods because it takes into consideration the dynamic 6 of 14 nature of the binding. In practice, the binding energy values obtained with both protocols are not comparable because they are obtained from different scoring functions. Therefore, itofwas the opted to select binding. the best In practice, the candidates in each binding energy protocol values for the obtained withfollowing steps.are both protocols So,not con- sidering comparablebothbecause strategies, theythe areproposed obtained hydrogel candidates from different scoringwere formed functions. by the following Therefore, it was blocks: opted to MN3-CM7-CL1 (protocol 1)inand select the best candidates eachMN4-CM5-CL1 protocol for the(protocol following2). For the steps. So, next analysis, considering both strategies, these hydrogels thewillproposed be named hydrogel candidates HG1 (derived fromwere formed protocol 1)by andtheHG2 following blocks: (derived from MN3-CM7-CL1 protocol 2). (protocol 1) and MN4-CM5-CL1 (protocol 2). For the next analysis, these hydrogels will be named HG1 (derived from protocol 1) and HG2 (derived from protocol 2). Figure 2. Ranking of MN. The boxplots of the cAMP-MNs binding energies calculated from the Figure 2. Ranking of MN. The boxplots of the cAMP-MNs binding energies calculated from the VINA scoring function are shown in (A,B) for the peptides MP-L and NCP-3a, respectively. The same VINA scoring function are shown in (A) and (B) for the peptides MP-L and NCP-3a, respectively. as (A,B) but for the binding energy obtained from the glide scoring function (IFd) are shown in (C,D). The same as (A) and (B) but for the binding energy obtained from the glide scoring function (IFd) **** correspond to p < 0.0001 of comparative statistical difference (using the t-test) with the lowest are shown in (C) and (D). **** correspond to p < 0.0001 of comparative statistical difference (using energy the building t test) block. with the lowest energy building block. 3.2. Molecular Dynamics Simulations In order to study the behavior of hydrogels and their interactions with cAMPs over time and in an environment that mimics reality, MDs were performed. Four systems were built, which were the combination of the two peptides MP-L (P1) and NCP-3a (P2) and the two proposed hydrogels (HG1 and HG2). It is worth mentioning that for a better understanding, the systems were labeled as follows: HG1_P1, HG1_P2, HG2_P1, and HG2_P2. First, the structural stability of the hydrogels and cAMP in the systems proposed was analyzed. Figure 5A shows the radius of gyration of the hydrogels during the 50 ns of simulation. Here, it can be observed that both HGs have similar behavior and that the radius of gyration does not vary too much between them, reaching values close to 40 Å. However, statistically significant differences are observed (Figure 5B). Both HGs have a larger radius of gyration in the presence of P1, and a smaller radius of gyration when in the presence of P2. This allows us to infer that small structural variations of the HGs occur to improve the interactions with each type of cAMP. Regarding the stability of the cAMP, Figure 6 shows the secondary structure elements (SSE) of each one of the cAMP’s residues in the four proposed systems. It is worth mentioning that each system (simulated Pharmaceutics 2023, 15, 474 7 of 14 in triplicate) was composed of 10 cAMPs, so each graph shows the statistical mode of the secondary structure for each residue at each specific time point. Overall, the results obtained show that the cAMPs maintain their secondary structure during the simulation. Figure 6A,B shows that the cAMP P1 retains its alpha-helical structure, and only small changes are observed between residues 12 to 15. In the case of the cAMP P2 (Figure 6C,D), it occurs in a similar way, where both beta sheets and the turn region in the center of the structure are maintained. From this last result, it can be inferred that both cAMPs7 do Pharmaceutics 2023, 14, x FOR PEER REVIEW of 14 not lose their secondary structure when interacting with the HGs, which is an important advantage for the cAMPs’ stability to maintain their antimicrobial function. Figure 3. Figure Ranking of 3. Ranking of CM. CM. The Theboxplots boxplotsofofthethecAMP-CMs cAMP-CMs binding bindingenergies calculated energies fromfrom calculated the VINA the scoring function are shown in (A,B) for the peptides MP-L and NCP-3a, VINA scoring function are shown in (A,B) for the peptides MP-L and NCP-3a, respectively. respectively. The same as The (A,B)as same but for the (A,B) butbinding energy obtained for the binding from the from energy obtained glide the scoring glidefunction scoring(IFd) are shown function in (C,D). (IFd) are shown in (C,D). * correspond * correspond to** to p < 0.05, p