Drug Repurposing by Simulating Flow Through Protein-Protein Interaction Networks.

As drug development is extremely expensive, the identification of novel indications for in-market drugs is financially attractive. Multiple algorithms are used to support such drug repurposing, but highly reliable methods combining simulation of intracellular networks and machine learning are currently not available. We developed an algorithm that simulates drug effects on the flow of information through protein-protein interaction networks, and used support vector machine to identify potentially effective drugs in our model disease, psoriasis. Using this method, we screened about 1,500 marketed and investigational substances, identified 51 drugs that were potentially effective, and selected three of them for experimental confirmation. All drugs inhibited tumor necrosis factor alpha-induced nuclear factor kappa B activity in vitro, suggesting they might be effective for treating psoriasis in humans. Additionally, these drugs significantly inhibited imiquimod-induced ear thickening and inflammation in the mouse model of the disease. All results suggest high prediction performance for the algorithm.

WHAT IS THE CURRENT KNOWLEDGE ON THE TOPIC?

☑ Several methods have been developed for drug repurposing, although these are mainly based on drug–disease associations and do not consider intracellular interactions. Additionally, only a few of them use machine learning for data processing.

WHAT QUESTION DID THIS STUDY ADDRESS?

☑ If it is possible to create a reliable drug repurposing algorithm with modeling intracellular protein–protein networks.

WHAT THIS STUDY ADDS TO OUR KNOWLEDGE

☑ With the use of intracellular network modeling and machine learning, it is possible to repurpose drugs with high reliability. Additionally, new interesting drug candidates and drug targets were identified in our study for the treatment of psoriasis.

HOW THIS MIGHT CHANGE CLINICAL PHARMACOLOGY OR TRANSLATIONAL SCIENCE

☑ With the increasing cost of drug development, drug repurposing is becoming a promising alternative. Our algorithm has the potential to significantly decrease the time and cost of drug candidate selection for a new indication.

Drug development is a long and extremely expensive process.1 The use of marketed drugs for new indications, referred to as drug repurposing, can significantly decrease expenses.1 Algorithms used for drug repurposing analyze several kinds of information either individually or in combination, including gene expression data, chemical properties of substances, drug target data, and static properties of intracellular networks.2, 3 Networks are reliable frameworks for the simulation of biological processes. Simulation of information flow through networks has already been used in drug discovery, but has yet to be used in the field of drug repurposing.4, 5, 6 We assumed that the use of such simulated dynamic data could potentially improve repurposing algorithms. As large datasets are generated during the simulation process, optimized dataset processing could also improve prediction performance. Machine learning, which is already widely used in drug development, is a clear choice for these optimization objectives. Different learning methods are available for predicting drug toxicity, potential inhibitors or activators of a given protein in compound libraries, pharmacokinetics, side effects, and drug interactions.7, 8, 9, 10 Support vector machine (SVM), one of the most widely used learning algorithms, classifies samples based on a model constructed from a training dataset. The model is used to differentiate between two groups of observations (two‐class SVM) or to identify observations similar to training data (one‐class SVM).11 We expected that the combination of network flow simulation and drug classification with SVM can be very effective in drug repurposing.

In our study we constructed intracellular protein–protein interaction networks, simulated information flow within these networks using a novel algorithm, and used SVM to identify potentially effective drugs in psoriasis. We also determined the in vitro and in vivo efficacy of identified drug candidates.

RESULTS

Network construction and flow simulation

A directed protein‐protein interaction (PPI) network was constructed from known interactions (activation or inhibition) between human proteins. A weighted adjacency matrix was generated and was used for the stepwise calculation of protein activity values. As basal activity values were added to the actual activity of all nodes before each step, the total activity of the system increased from step to step when all interactions in the network were activating. Because there were sufficient inhibiting interactions in the constructed network, equilibrium (defined as a change of <1% in the total activity of the system between two subsequent steps) was achieved.

Initially, the algorithm was executed using unit basal activity of all nodes until equilibrium was achieved. Using data from our recent microarray meta‐analysis,12 basal activity values at equilibrium were then adjusted to reflect conditions in psoriatic lesional skin: basal activities of protein‐coding genes found to be upregulated in lesional skin samples were doubled and the activity of protein‐coding genes found to be downregulated were halved. Simulation was executed with the new basal activity values until a new equilibrium was achieved.

Simulation of drug effects and prediction of drug efficacy

Drug effects were simulated with different hypothetical drug efficiency values (ed in Eqs. 2 and 3). We also simulated drug effects by changing the basal activity of target nodes only, without changing the flow capacity of the corresponding interactions. For each drug, basal activity values of target nodes and flow capacity of the corresponding interactions were applied to the equilibrium values achieved for psoriatic skin conditions and drug effects were simulated until equilibrium was achieved. Drug‐ and efficiency‐specific spreading matrices were obtained for all 1,440 drugs simulated (Figure 1). SVM training was performed with actual node activity values at equilibrium for 18 drugs indicated in the treatment of psoriasis (Supplementary Table 1). Training was performed with and without dimension reduction of data and with different kernels (Table 1). Only models with higher than 50% mean cross‐validation accuracy were used for effective drug prediction. To avoid false‐positive results, drugs predicted to be effective with all reliable models were considered.

Figure 1

Method for predictive repurposing of drugs for a given indication. Network information flow is simulated in intracellular networks and drug‐specific spreading matrices (i.e., matrices containing protein activity values for different steps) are generated. SVM is trained with protein activity data for drugs already used for the given indication. SVM models are used to predict effective drugs. False‐positive results are eliminated by multiple simulation and model construction steps. Drugs that are predicted to be effective in most high‐accuracy SVM models are considered potentially effective. [Color figure can be viewed at cpt-journal.com]

Table 1

Properties for the training of SVM

Drug efficacy	Dimension reduction	Kernel	Gamma	μ	Degree	Accuracy	Used for prediction
100% (only node)	No	Linear	NA	12.1	NA	85.7143	Yes
100% (only node)	No	Polynomial	3	–18.1	3	64.28575	Yes
100% (only node)	Yes	Linear	NA	–8.1	NA	64.28575	Yes
100% (only node)	Yes	Polynomial	–1	–18.1	3	57.1429	Yes
100% (only node)	Yes	Polynomial	–15	–0.1	10	35.7143	No
100%	No	Linear	NA	–12.1	NA	85.7143	Yes
100%	No	Polynomial	–1	–18.1	3	71.4286	Yes
100%	Yes	Linear	NA	–8.1	NA	71.4286	Yes
100%	Yes	Polynomial	3	–16.1	3	71.4286	Yes
100%	Yes	Polynomial	–1	–20.1	10	42.85715	No
50%	No	Linear	NA	–10.1	NA	85.7143	Yes
50%	No	Polynomial	–1	–18.1	3	64.28575	Yes
50%	Yes	Linear	NA	–14.1	NA	71.4286	Yes
50%	Yes	Polynomial	3	–10.1	3	71.4286	Yes
50%	Yes	Polynomial	–1	–20.1	10	50	No
20%	No	Linear	NA	–18.1	NA	85.7143	Yes
20%	No	Polynomial	–1	–20.1	3	64.28575	Yes
20%	Yes	Linear	NA	–12.1	NA	71.4286	Yes
20%	Yes	Polynomial	3	–18.1	3	64.28575	Yes
20%	Yes	Polynomial	–15	–0.1	10	57.14285	Yes

Different drug efficacies were simulated. Training of SVM was carried out with dimension‐reduced data as well and different kernel types. Third and tenth degree polynomials were used for polynomial kernels. Gamma and μ parameters were assessed by parameter optimization (gamma is not applicable for linear kernel). Models with the highest accuracy generated during parameter optimization (at least 50%) were used for prediction.

Fifty‐one drugs were predicted to be effective by all high‐accuracy models (Table 2). Based on drug actions, the five most prevalent drug groups were gamma‐aminobutyric acid (GABA) receptor activators (barbiturates), sympathomimetic drugs, vitamin D receptor (VDR) agonists, GnRH‐receptor agonists, and PPARG agonists. Barbiturates are yet not associated with psoriasis, but immunomodulatory effects of these substances have been suggested.13 Epinephrine is presumed to have a role in psoriasis pathogenesis, and, thus, the modulation of adrenergic signaling may have disease‐modifying effects.14, 15 As expected, VDR agonists, which are first‐line therapeutic options for psoriasis, were predicted to be effective.16 GnRH is not yet associated with psoriasis and the role of GnRH signaling has not yet been investigated for this disease. Although PPARG agonist thiazolidinediones have not yet been indicated for the treatment of psoriasis, several clinical studies prove beneficial disease‐modifying effects of these drugs.17 We chose three drugs that previously have not been directly associated with psoriasis and for which no in vitro and/or in vivo data are available regarding their disease‐modifying effects to confirm our prediction results: the barbiturate hexobarbital sodium, the β2‐adrenergic agonist salbutamol hemisulphate, and the GnRH‐receptor agonist leuprolide acetate (LPA).

Table 2

Potentially effective drugs for the treatment of psoriasis

Classification	Count	Drugs
GABA(A) activator (theoretical data)	17
N05CA Barbiturates, plain	11	Amobarbital, Aprobarbital, Barbital, Barbituric acid derivative, Butabarbital, Heptabarbital, Hexobarbital, Pentobarbital, Secobarbital, Talbutal, Thiopental
N03AA Barbiturates and derivatives	4	Metharbital, Methylphenobarbital, Phenobarbital, Primidone
Nonclassified	2	Butalbital, Butethal
Symphatomimetic (theoretical data)	15
R03AC Selective beta‐2‐adrenoreceptor agonists	6	Indacaterol, Pirbuterol, Procaterol, Salbutamol, Salmeterol, Terbutaline
R03CC Selective beta‐2‐adrenoreceptor agonists	3	Bambuterol, Fenoterol, Formoterol
C01CA Adrenergic and dopaminergic agents	1	Arbutamine
R03AB Nonselective beta‐adrenoreceptor agonists	1	Orciprenaline
R03CA Alpha‐ and beta‐adrenoreceptor agonists	1	Ephedra
R03CB Nonselective beta‐adrenoreceptor agonists	1	Isoproterenol
G02CA Sympathomimetics, labor repressants	1	Ritodrine
S01EA Sympathomimetics in glaucoma therapy	1	Dipivefrin
Nonclassified	1	Arformoterol
VDR agonist (studies available)	6
A11CC Vitamin D and analogues	5	Alfacalcidol, Calcidiol, Cholecalciferol, Dihydrotachysterol, Ergocalciferol
H05B Anti‐parathyroid agents	1	Paricalcitol
GnRH Agonist (not investigated)	4
H01CA Gonadotropin‐releasing hormones	2	Gonadorelin, Nafarelin
L02AE Gonadotropin releasing hormone analogues	2	Goserelin, Leuprolide
PPARG agonist (studies available)	2
A10BG Thiazolidinediones	2	Pioglitazone, Rosiglitazone
Adenosine receptor agonist (clinical trials)	1
C01E Other cardiac preparations	1	Regadenoson
COX inhibitor (controversial)	1
M01A Anti inflammatory and antirheumatic products, nonsteroids	1	Ibuprofen
Potassium channel opener (not investigated)	1
C02D Arteriolar smooth muscle, agents acting on	1	Minoxidil
Dopaminergic agonist (pilot study)	1
C01CA Adrenergic and dopaminergic agents	1	Dobutamine
Nonproteinogenic amino acid (not investigated)	1	Canaline
Nonsaccharide sweetener (not investigated)	1	Aspartame

Drugs were classified based on their action and association of the drug class with psoriasis is indicated in parentheses. Subclassification was based on the Anatomical Therapeutic Chemical (ATC) Classification System.

In vitro confirmation of potentially effective drugs

To select the ideal in vitro setup for the examination of drug actions, we carried out analysis of the model PPI network to identify the most relevant psoriasis‐associated cytokines and pathways, which are related to the action of predicted drugs. We calculated the betweenness node centrality value of proteins in the model PPI network and in a mixed network containing both protein–protein interactions and interactions between predicted drugs and their targets. Betweenness is the number of shortest paths crossing through a given node in the network. We assumed that most relevant nodes related to drug action are the ones that betweenness centrality increased the most in the mixed network compared to the model PPI network. Tumor necrosis factor alpha (TNF) was the top‐ranked cytokine and nuclear factor kappa B (NFκB) was above the 98th percentile in the results (Supplementary Table 2). This cytokine and activated signaling pathway plays an essential role in the pathogenesis of psoriasis. Additionally, most drugs indicated for the treatment of psoriasis alter NFκB activity and the examination of this pathway is the focus of drug discovery for the disease.18, 19, 20, 21 Thus, the effect of different drug concentrations on TNF‐induced NFκB activity was measured in our experiments using a dual‐luciferase reporter assay. As hyperproliferation of keratinocytes is one of the most important hallmarks in psoriasis, the TNF‐induced NFκB activity was measured in a keratinocyte cell line.22 A plasmid containing NFκB promoter‐controlled firefly luciferase gene and a plasmid containing a renilla luciferase gene were both transfected into HPV‐keratinocyte cells. Cells were treated with three different concentrations of hexobarbital‐sodium and salbutamol‐hemisulphate 24 h after transfection, and 10 ng/ml human recombinant TNF was added to the cultures 1 h after treatment. Relative luminescence, which is proportional to NFκB induction, was measured 17 h after TNF treatment. This experimental setup produced inconsistent results for leuprolide‐acetate, which may be explained by the finding that GnRH effects are mediated by the long‐term pulsatile release of the hormone.23 Thus, three different concentrations of leuprolide‐acetate were used for 72 h, replacing the drug‐containing media every 24 h. We added 10 ng/ml human recombinant TNF after 72 h and measured luminescence 17 h later.

TNF‐induced NFκB activity was suppressed by all the tested drugs, where inhibition was dose‐dependent for hexobarbital‐sodium and leuprolide‐acetate (Figure 2). To further examine whether suppression of the NFκB activity results in modulation of the pathway, we carried out quantitative real‐time polymerase chain reaction (qRT PCR) for the NFκB target gene interleukin (IL)‐1b. All tested drugs inhibited the upregulation of IL‐1b after TNF treatment (Supplementary Figure 1). Our results suggested potential disease‐modifying effects of these drugs in vitro. Considering the role of keratinocytes in psoriasis, we next examined the effects of these drugs on lesion formation in a mouse model of psoriasis.

Figure 2

All drugs predicted for efficacy in psoriasis significantly inhibited the TNF‐dependent induction of NFκB in HPV‐keratinocytes. Data are presented as means ± SD. *P < 0.001 vs. control group; **P < 0.001 vs. 10 ng/ml TNF group; ***P < 0.01 vs. 10 ng/ml TNF group; n = 9. HXB, hexobarbital sodium; SBT, salbutamol‐hemisulphate; LPA, leuprolide‐acetate. [Color figure can be viewed at cpt-journal.com]

In vivo confirmation of predicted drugs

We carried out in vivo confirmation of potentially effective drugs. Imiquimod‐induced skin inflammation in mice is a conventional model of psoriasiform dermatitis.10 We treated mice for 6 days with imiquimod (positive control), imiquimod and gels containing drugs, or imiquimod and gel vehicle alone (negative control). Ear thickness, which is indicative of edema, cellular infiltration, and hyperproliferation of the epidermis, all of which are characteristic of psoriasis, was measured each day. All three drugs significantly reduced ear thickening (Figures 3, 4 a): whereas ear thickness increased 80–90% in mice treated with imiquimod and imiquimod + vehicle compared to untreated mice, it increased only 20–30% in drug‐treated mice. Histological examination of ear specimens was performed after 6 days of treatment. As keratinocyte hyperproliferation is characteristic for psoriasis, we measured epidermal thickness of all ear specimens at multiple loci, and found that ear thickness was significantly lower in drug‐treated mice compared to positive and negative control mice (Figure 4 b). We also assessed cellular infiltration in the dermis computationally. Cellular infiltration is caused by skin homing of T‐lymphocytes and extravasation of neutrophil granulocytes induced by inflammatory cytokines.16 We found that all drugs significantly decreased cellular infiltration by Day 6 compared to positive and negative control treatment (Figure 4 c).

Figure 3

All predicted drugs significantly inhibited imiquimod‐induced ear thickening in mice. (a) Ear thickness (μm). (b) Relative ear thickness. All measured values are normalized to the mean ear thickness for the group on day 0. Data are presented as means ± SD. *P < 0.05 vs. measurements for mice treated with imiquimod or imiquimod + vehicle; n ≥ 10. IMQ, imiquimod; HXB, hexobarbital sodium; SBT, salbutamol‐hemisulphate; LPA, leuprolide‐acetate.

Figure 4

Histological examination of ear specimens after 6 days of treatment. (a) Characteristic microscopic pictures from each treatment group, (b) thickness of the ear epidermis, and (c) average cell count. All three datasets indicated that the drugs significantly decreased the thickening of epidermis and cellular infiltration. Data are presented as mean ± SD; *P < 0.05; n ≥ 10. IMQ, imiquimod; HXB, hexobarbital sodium; SBT, salbutamol‐hemisulphate; LPA, leuprolide‐acetate.

DISCUSSION

Analysis and modeling of intracellular networks are popular methods in drug discovery, and our results suggest that these methods, in combination with machine learning, are powerful tools for drug repositioning as well. We developed a promising algorithm to simulate drug effects on protein–protein interaction networks. We selected psoriasis to test our method, as this common disease currently has multiple treatment options. We generated multiple high‐accuracy models with SVM for drug efficacy and used these models to predict potential disease‐modifying substances. To avoid false‐positive results, we applied stringent selection criteria to the prediction results.

To validate the algorithm we selected three drugs with potential antipsoriatic activity. Barbiturates activate the GABA receptors and are used as sedatives and anticonvulsants. Historically, skin symptoms of barbiturate poisoning already suggested an association of these substances with the skin.24 GABA type “A” receptors were identified in the epidermis of mice and activation of this receptor prevented epidermal hyperplasia.25 The effect of GABA signaling on immune functions has already been investigated. Thiopental was reported to inhibit the activity of IκB kinase and, thus, NFκB in Jurkat T‐lymphocytes. However, this effect was attributed to the thio‐group and was suggested to be specific for this substance.26 Our prediction and experimental validation suggested a positive effect of hexobarbital on psoriasiform skin inflammation. We assume that the inhibition of NFκB may be mediated by the increased intracellular calcium level caused by GABA receptor activation, as calcium is already reported to inhibit the NFκB pathway in some cell types.27

Salbutamol is a β2 adrenergic receptor agonist indicated for the treatment of bronchospasm inter alia.28 The presence of such receptors on the membrane of keratinocytes has already been demonstrated.29 The polymorphism of such receptors is associated with susceptibility to psoriasis, which suggests a potential effect of its modulation on disease activity.30 The β2 adrenergic signaling pathway is associated with several immune‐related functions and with the modulation of NFκB. However, the mode of action and whether pro‐ or antiinflammatory mediators are secreted are highly dependent on cell type. The alteration of the NFκB pathway by beta adrenergic signaling is reported to be expressed through several mechanisms, like the modulation of NFκB subunits or interaction between membrane receptors of the two signaling pathways.28 Our results indicated salbutamol‐mediated NFkappaB inhibition in vitro and its antipsoriatic effect in vivo.

Leuprolide acetate activates the GnRH receptor and is indicated for the treatment of prostate cancer and endometriosis, among others.31 The presence of the GnRH receptor in human keratinocytes has not been reported yet, but immunohistochemical studies of canine skin reported its presence in the epidermis.32 Several reports are available about the effects of GnRH on the immune system. Modulation of NFκB was also reported by GnRH in a macrophage cell line.33 In that study, GnRH treatment resulted in NFκB inhibition, which may be caused by increased intracellular calcium levels. Although a case report described precipitation of psoriasiform eruptions during LPA therapy, our results suggest a potentially beneficial effect of this substance on psoriasis.34

Although additional clinical evaluation is needed to prove the efficacy of these drugs on psoriasis in humans, the security assessment for these substances has been completed in previous clinical studies, potentially reducing the costs of repurposing these drugs for psoriasis. Our models were constructed using simulation results for drugs that are already used in the treatment of psoriasis. Although we obtained promising results in psoriasis, in the case of rare diseases or diseases with fewer treatment options we would expect our method to produce models with lower reliability. Our in vitro and in vivo results suggested that the selected drugs would be effective in treating psoriasis, but our method did not predict the level of drug efficacy. Since the prediction is based on protein–protein interactions, the use of the method is limited to drugs that mediate protein activities. Future testing and optimization of our method should include other diseases and drugs to determine its general applicability.

MATERIALS AND METHODS

Network flow simulation

We developed an algorithm in R (Vienna, Austria) that calculates stepwise activity of nodes in a network. One step in the flow is the diffusion of all actual activity from one node (source) to downstream connected nodes (targets). The model uses a weighted adjacency matrix that is generated from the edge list format of networks. All nodes have a basal activity added to its actual activity before the next step. In a step, the activity flows from the source node “a” to the target node “b”, as follows: where wa,b = the weight of the interaction (positive real number), ma,b = action of node “a” on “b” (m = 1 in case of activation, m = –1 in case of inhibition), xa = activity of the source node before actual step, c = basal activity of source node, and n = number of interactions originating from node “a”.

The algorithm creates a so‐called spreading matrix, in which rows correspond to nodes and columns correspond to actual activity values of nodes at the beginning of a given step. The spreading matrix is created using adjacency matrix and Eq. 1.

Machine learning and simulation of drug effects

Drug effects were simulated by changing the basal activity of drug target proteins and/or the flow capacity through interactions belonging to the target node. The effect of 1,440 in‐market or experimental drugs was simulated and drug‐specific spreading matrices were generated. Effects of a drug “d” on flow through interactions of target node “a” was calculated with the following formulas: where a = drug target, b = target node of node “a” f = node acting on node “a” e = drug efficacy (real values between 0 and 1), g = mode of drug action on target node (1 for activation; –1 for inhibition), w = the weight of the interaction (positive real number), m = action type (m = 1 for activation, m = –1 for inhibition), x = activity of the source node before actual step, c = basal activity of source node, n = number of interactions originating from source node.

The e1071 R library was used for machine learning. This package offers an interface to the LIBSVM library.35, 36 We selected one‐class SVM, which differentiates one particular class of observations from all other observations (for example, effective drugs from all other drugs). Feature scaling and parameter optimization was carried out before SVM training. Parameter optimization was performed with 2‐fold cross‐validation. Multiple SVM runs were carried out with linear, third, and tenth degree polynomial kernels. All runs were repeated with dimension reduction using principal component analysis in R. Dimensions were reduced so that the loss of variance was less than 1%.

Online data sources

The list of protein–protein interactions was downloaded from the STRING 9.0 database.37 Highly reliable (confidence score greater than 0.800) directed interactions available for Homo sapiens with known action (i.e., activation or inhibition) were used for protein–protein interaction network construction. Drug–protein interaction data with known mode of action (activation or inhibition) were downloaded from Drugbank.31 Drugs indicated for the treatment of psoriasis were collected from the US Food and Drug Administration (FDA) Drug Label Database (https://rm2.scinet.fda.gov/druglabel#simsearch‐0). KEGG DRUG was used for the ATC classification of drugs.38

Cell culturing

The HPV‐keratinocyte cell line (HPV‐KER, established in our laboratory from normal human keratinocytes by Polyanka et al.) was cultured in serum‐free complete keratinocyte medium (KSFM, Gibco Laboratories, Grand Island, NY) supplemented with 5 ng/ml recombinant epidermal growth factor, 50 μg/ml bovine pituitary extract, L‐glutamine, and antibiotic/antimycotic solution in a CO2 atmosphere at 37°C. Experiments were carried out with cells that had undergone at least two passages after thawing.

Dual luciferase reporter assay

HPV‐keratinocytes (5 × 105) were plated in 12‐well plates from flasks containing cells at 70–80% confluence. Culturing media was refreshed before transfection. Cells at 70% confluence were transfected with the pGL4.75[hRluc/CMV] vector, which contains the hRluc reporter gene regulated by a cytomegalovirus (CMV) immediate‐early enhancer/promoter, and the PathDetect NF‐κB cis‐Reporting plasmid (pNF‐κB‐Luc, Agilent Technologies, Palo Alto, CA), which contains direct repeats of the transcription recognition sequences for NF‐κB, that drives transcription of the luc2P luciferase reporter gene. The transfection mix contained 25 ng/ml pGL4.75, 10 μg/ml pNF‐κB‐Luc, 3 V/V% Roche, Basel, Switzerland, x‐tremeGENE9 DNA transfection reagent in 100 μl KSFM and was added to wells containing 500 μl media. Twenty‐four hours after transfection, media were removed and cells were first washed with phosphate‐buffered saline (PBS, Lonza, Basel, Switzerland) and then were treated with drug candidates: hexobarbital‐sodium (VEB Arzneimittelwerk, Dresden, Germany), salbutamol‐hemisulphate (Sigma‐Aldrich, St. Louis, MO), and leuprolide‐acetate (Sigma‐Aldrich). Drugs were dissolved in KSFM to three different concentrations. One hour after treatment, 10 ng/ml TNF (R&D Biosystems, Minneapolis, MN) was added to the media. For chronic leuprolide‐acetate treatment, LPA‐containing media was used for 72 h and refreshed every 24 h. Ten ng/ml TNF was added to the media after 72 h of treatment. Culturing media was removed 17 h after TNF‐treatment, cells were washed with PBS and 500 μl lysis buffer (100 μl of Promega, Madison, WI, Reporter 5× lysis buffer in 400 μl of distilled water) was added to each well. Plates were shaken for 15 min and the contents of each well were collected and centrifuged at 10,000g for 30 sec. Supernatants were used for assay. All experiments were carried out in three technical and three biological replicates; each biological replicate was performed with new, freshly thawed HPV‐keratinocyte pools. Luminescence indicating firefly and renilla luciferase activity was measured with a Thermo Scientific (Pittsburgh, PA) Luminoskan Ascent Microplate Luminometer. All samples were measured three times and average values were used for further calculation. Background values of plate wells were measured and subtracted from sample values. Relative luciferase unit values for samples were calculated by dividing firefly relative luminescence unit values with renilla relative luminescence unit values. A Shapiro–Wilk normality test was used to check whether measured values follow normal distribution. The identification of relevant pathways for in vitro experiments was carried out with network analysis using the igraph “R” package.39

RNA isolation and qRT PCR

For RNA isolation, cells were suspended in 500 μl TRIzol reagent (Sigma‐Aldrich) and incubated for 10 min. Then 100 μl chloroform was added and the samples were vortexed. The samples were centrifuged for 20 min at 10,000g and the upper aqueous phase was transferred to a clean tube. The RNA was precipitated by adding 250 μl isopropanol (Sigma‐Aldrich). The samples were incubated for 10 min and then centrifuged for 20 min at 10,000g. The supernatant was discarded and the RNA pellet was resuspended in 75% ethanol. After 15 min incubation, the samples were centrifuged for 5 min at 10,000g. The ethanol was discarded and the RNA pellet was dried at room temperature. The RNA pellet was dissolved in 20 μl DEPC‐treated water. cDNA was synthesized from 0.5 μg RNA using RevertAid First Strand cDNA Synthesis Kit (Thermo Scientific) following the manufacturer's protocol. Abundance of IL‐1b and 18s RNA transcripts was quantified using custom primer sets and the Universal Probe Library (Roche) with qPCRBIO Probe Mix (PCR Biosystems, London, UK) reagent. Relative gene expression was calculated by first normalizing data for the 18S ribosomal RNA expression and then using the 2‐ΔΔCt method. In the case of luciferase assay and qRT PCR results, analysis of variance (ANOVA) was used for the comparison of group variances; a pairwise t‐test with Bonferroni correction was used to compare groups. Between‐group difference with corrected P < 0.05 was considered significant. The following primers were used for qRT PCR: 18s 5′: CGCTCCACCAACTAAGAACG; 18s 3′: CTCAACACGGGAAACCTCAC; IL1b 5′: AAAGCTTGGTGATGTCTGGTC; IL1b 3′: AAAGGACATGGAGAACACCACT.

Preparation of gel

Active pharmaceutical ingredient‐free hydrogels were prepared with 1.56 wt% hydroxypropyl methylcellulose (METHOCEL E4M Premium, Dow Chemical, Midland, MI) and distilled water. Hexobarbital sodium (1.3 wt%), salbutamol hemisulfate (0.42 wt%), or leuprolide acetate (0.38 wt%) was dissolved in the hydrogel formulation. The pH of the final hydrogels were 6.0, except the hydrogel containing LPA, which was adjusted to a pH of 5.2 by adding an aqueous citric acid solution (Hungaropharma, Hungary) to ensure the stability of the LPA.40

In vivo experiments

Experiments were carried out on 8–9‐week‐old male BALB/c mice weighing 21–25 g. The animals were housed in plastic cages in a thermoneutral environment with a 12‐h light‐dark cycle and had access to standard laboratory chow and water ad libitum. All interventions were in full accordance with the National Institutes of Health (NIH) guidelines and protocols were approved by the Ethical Committee for the Protection of Animals in Scientific Research at the University of Szeged (license number: I.74‐3/2015.MÁB).

On the day of the experiment, initial thickness of both ears of all mice was measured with a thickness gauge (“Oditest,” Kroeplin, Germany). Ear thickness measurement was carried out on all subsequent days before any treatment. All treatment groups were treated with 5% imiquimod cream (Aldara–MEDA Pharma, Hungary) at 13:00 every day. All groups except the positive control group were treated with drug or vehicle at 9:00 and 17:00 every day. After 6 days of treatment, mice were euthanized with an overdose of ketamine (300 mg/kg) and histological samples were taken from both ears.

Histology

Ears of mice were fixed in a buffered solution of formaldehyde (4%) and embedded in paraffin from which 3‐μm‐thick coated slides were prepared. All slides were stained with hematoxylin and eosin. Microscopic images of histological specimens were digitalized and visualized with Pannoramic Viewer v. 1.15. (3DHISTECH). Epidermal thickness of all ear specimens was measured 10 times in all field of views at 10× magnification by two independent persons and the mean thickness was calculated for each ear specimen. The number of cells in the dermis was calculated computationally using the ImageJ software41 in all homogenous fields of view (i.e., excluding vessels, sebaceous glands, epidermal, subcutaneous structures, and follicles) at 40× magnification for each ear specimen. Briefly, the areas to analyze were converted to 8‐bit color‐depth images and the color was changed to black and white based on a threshold. Subsequently, watershed segmentation was carried out, and particles were analyzed with a size property of 500 to infinity and circularity property 0 to 1. Mean cell count in a field of view was calculated for each ear specimen.

CONFLICT OF INTEREST

US nonprovisional patent application is in process for the algorithm (M. Manczinger, L. Kemény ‐ US 15/198,366).

AUTHOR CONTRIBUTIONS

M.M. and L.K. wrote the article; M.M., G.E., and L.K. designed the research; M.M., V.B., S.B.B., K.S., E.S., and G.E. performed the research; M.M., V.B., and B.T.P. analyzed the data; M.M., B.B., and E.C. contributed new reagents/analytical tools.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.

Supporting Information

Click here for additional data file.