• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    PIMD:An Integrative Approach for DrugRepositioning Using Multiple Characterization Fusion

    2020-09-02 03:59:14SongHeYuqiWenXiaoxiYangZhenLiuXinyuSongXinHuangXiaochenBo
    Genomics,Proteomics & Bioinformatics 2020年5期

    Song He,Yuqi Wen,Xiaoxi Yang,Zhen Liu,Xinyu Song,Xin Huang,Xiaochen Bo*

    Department of Biotechnology,Beijing Institute of Radiation Medicine,Beijing 100850,China

    KEYWORDS Drug repositioning;Drug similarity network;Multiple characterization fusion;Network pharmacology;Drug discovery

    Abstract The accumulation of various types of drug informatics data and computational approaches for drug repositioning can accelerate pharmaceutical research and development.However,the integration of multi-dimensional drug data for precision repositioning remains a pressing challenge.Here,we propose a systematic framework named PIMD to predict drug therapeutic properties by integrating multi-dimensional data for drug repositioning.In PIMD,drug similarity networks(DSNs)based on chemical,pharmacological,and clinical data are fused into an integrated DSN(iDSN)composed of many clusters.Rather than simple fusion,PIMD offers a systematic way to annotate clusters.Unexpected drugs within clusters and drug pairs with a high iDSN similarity score are therefore identified to predict novel therapeutic uses.PIMD provides new insights into the universality,individuality,and complementarity of different drug properties by evaluating the contribution of each property data.To test the performance of PIMD,we use chemical,pharmacological,and clinical properties to generate an iDSN.Analyses of the contributions of each drug property indicate that this iDSN was driven by all data types and performs better than other DSNs.Within the top 20 recommended drug pairs,7 drugs have been reported to be repurposed.The source code for PIMD is available at https://github.com/Sepstar/PIMD/.

    Introduction

    Despite the ever-increasing funding of pharmaceutical research and development (R& D),the number of new drugs approved has not increased significantly [1].Traditionalde novodrug development remains costly,risky,and time-consuming [2,3].Drug repositioning,wherein an existing drug receives a new application,provides a new opportunity for pharmaceutical R& D[4].With the accumulation of drug informatics datasets,computational algorithms can be used for systematic identification of potential new indications for on-market drugs,thereby reducing the financial and time investment,as well as the risk,involved in pharmaceutical R& D [5-8].

    Recent studies have shown that computational approaches based on drug similarity have the potential to reveal novel indications for on-market drugs.These approaches have been applied primarily in the following four ways.1)The first group includes transcriptional response-based approaches.For instance,Iorio et al.[9] used drug-specific response profiles based on the Connectivity Map(CMap)database to find drugs with similar modes of action,whereas Xie et al.[10]used drug perturbation profiles from the Library of Integrated Networkbased Cellular Signatures (LINCS) project to predict additional therapeutic properties of drugs.2) The second group includes chemical structure-based approaches.For example,Keiser et al.[11] predicted potential drug targets based on combined drug-target structure.3) The third group includes side effect-or other phenotype-based approaches.For example,Campillos et al.[12,13] constructed a side effect-driven drug similarity network (DSN) based on the assumption that drugs with similar side effects may share targets leading to the identification of novel drug targets.4) The last group includes target property-based approaches.For example,Yildirim et al.[14] used known drug-target associations to assess ongoing trends and shifts in drug discovery and to quantify interrelationships between drug targets and diseasecausing gene products.Although they shared the assumption that similar drugs tend to share therapeutic properties,each of these studies focused on a single drug feature in assessment of drug similarity,raising some doubts about the usefulness of these approaches.For example,Yildirim et al.[14]pointed out that most drugs with the same targets have different chemical structures,and Keiser et al.[15] demonstrated that a small change in drug structure could alter binding affinity dramatically.In addition,the transcriptional response to drug perturbation may differ across cell lines and drug dosages,thus introducing noise into drug repositioning strategies based on transcriptome data.Notwithstanding,our previous studies illustrated positive correlations between repositioning potential based on transcriptome data and that based on side effect profile or structure [10].

    With the development of network pharmacology and systems biology,integrating multi-attribute data of drugs seems to be a feasible means of identifying new opportunities for drug repositioning [16].One of the most common methods to integrate such data is to concatenate several measurements from various properties,such as side effects and chemical fragments,of each drug [17,18].However,the already low signalto-noise ratio in each data type could be diluted by concatenation [19].To avoid this problem,many researchers have made some preliminary attempts to use DSNs based on different drug properties for data combination.For example,Napolitano et al.[20] constructed three DSNs,based on drug structures,distances between drug targets in protein-protein interaction (PPI) networks,and expression patterns of drug perturbations,separately.They then integrated these attribute datasets by averaging three drug similarity measurements to predict new therapeutic properties of drugs.Meanwhile,Wang et al.[21]proposed a new algorithm,called PreDR,which predicts as yet unidentified drug-disease associations by taking the maximums of three drug similarity matrices derived from chemical structure,target protein sequence,and side effect profile similarities.Zhang et al.[22] proposed the Similaritybased LArge-margin learning of Multiple Sources (SLAMS)algorithm of drug similarity based on multiple sources of drug and disease property data.SLAMS outputs therapeutic scores for each drug-disease pair that correspond to multi-level drug properties and disease properties,and then averages the scores to predict the novel disease applications for drugs.Liu et al.[23] proposed the two-pass random walks with restart on a heterogeneous network (TP-NRWRH) to predict new indications for approved drugs.In the model,DSNs are integrated using the probability disjunction formula.Additionally,there are many articles on the drug-target predictions with data integration based on a linear combination of multiple attributes[24-42].Although many studies take into account multiple drug sources,these integration strategies based on simple averaging or maximization are linear and cannot make full use of topology and non-linear information.Furthermore,most of the methods cannot evaluate the contribution and the relative importance of each property data.As for network-based analysis,to our knowledge,there is no integrated method offering a systematic way to annotate and evaluate the drug clusters.In recent years,non-linear multi-dimensional data fusion algorithms and tools have been widely applied in disease subtype identification,such as the Similarity Network Fusion (SNF)and Integrated Clustering of Multidimensional biomedical data (ICM) [19,43].

    In this study,we propose a systematic and extensible paradigm of drug repositioning,namely prediction of drug therapeutic property by integrating multi-dimensional data(PIMD),and report the construction of an integrated DSN(iDSN) based on drug structure,side effect profile,and target protein sequence data.First,we integrated different types of drug information,including side effects,chemical structures,and molecular targets representing clinical,chemical,and pharmacological properties,separately,and constructed a DSN for each of these properties.Second,we used a nonlinear fusion algorithm,namely SNF,to combine three DSNs into the iDSN iteratively[19].Next,we evaluated the contributions of each dimension of data in the iDSN and the correlation among them.

    Our study examined the types of data underlying the iDSN and how the iDSN performs relative to single-property networks.We used spectral clustering to divide the iDSN into clusters,and then conducted a systematic and comprehensive drug cluster analysis through five types of statistical analyses,including drug-based enrichment analysis,target-based enrichment analysis,drug property analysis,chemogenomic enrichment analysis (CGEA),and chemical ontology enrichment analysis.We hypothesize that if similar drugs have similar therapeutic properties,then the drugs that appear unexpectedly in a cluster based on their anatomical therapeutic chemical (ATC) label would represent repositioning candidates.

    Method

    Data source

    We collected information on 7132 drugs with their corresponding target protein information from the DrugBank database,a bioinformatics resource with detailed drug data and complete drug-target interaction data.We used version 5.0 of the Drug-Bank database to construct the DSN based on drug targets(DSN-T) [5].Protein sequence data were extracted from the UniProt database,which provides high-quality,freely accessible protein sequence data [44].

    We obtained 139,756 relationships between 1430 drugs and 5868 side effects from the Side Effect Resource(SIDER)database,which contains information about on-market drugs and their recorded adverse drug reactions or side effects,extracted from public documents and package inserts.We used version 4.1 of SIDER to construct the DSN based on drug side effects(DSN-S) [8].

    PubChem Compound is a database containing more than 92 million unique structures of compounds.Similarly,We extracted chemical structures of drugs from PubChem Compound to construct the DSN based on drug chemical structure(DSN-C) [7].We used the PubChem Compound Identifier(CID) as the only identifier of drugs to identify drugs shared across the databases.We used the identified drugs to construct three single DSNs.

    Drug similarity measurements

    Drug similarity quantifies the degree of shared features between paired drugs.We restricted similarity scores to be between 0 (lowest) and 1 (highest).We defined the drug similarity measurements of the three properties examined separately.

    Drug similarity based on drug side effects

    Side effects represent the clinical properties of a drug.We obtained side effect data from SIDER.We used drug side effect information with frequency data(as opposed to without)because such information was derived empirically and thus deemed more credible.Given there are risks of bias in observation and statistics,we filtered outsider effect terms if they occurred only once or with a frequency <0.1%.Finally,we characterized drug side effects according to the 2072-dimensional binary vectorE(d),known as side effect profile.Similarity based on side effects between two drugsdandd′was computed by the Tanimoto coefficient of their side effect profiles:

    where |E(d)|and |E( d′)|are the number of side effect terms for drugs d and d′, respectively. E(d)×E(d′) represents the number of side effects shared by these two drugs.

    Drug similarity based on drug chemical structure

    Drug chemical structure represents the chemical properties of a drug.We obtained chemical structure information from Pub-Chem Compound and computed atom-pair descriptors of drugs using the R package ‘‘ChemmineR”[45].The atompair descriptors used to quantify the chemical structure of small molecule compounds encode all atom pairs in a drug.We computed similarity based on the chemical structure between drugsdandd′as the Tanimoto coefficient of the chemical atom-pair descriptors:

    where |C(d)|and |C( d′)|represent the number of atom pairs for drug d and drug d′, respectively. C(d)×C(d′) represents the number of the atom pairs shared by these two drugs.

    Drug similarity based on drug targets

    Drug targets represent the pharmacological properties of a drug.We obtained drug-target interactions and target protein sequences from the DrugBank and Uniprot databases,respectively.Then,similarity based on drug targets between drugsdandd′was computed with a normalized Smith-Waterman score as follows:

    where P(d)represents a target protein set of drug d,Pi(d)indicates the ith target of drug d, and |P(d)|represents the size of the target protein set of drug d.SW(Pi(d),Pi(d′))is the Smith-Waterman sequence alignment score of target proteins of drugs d and d′[46].

    SNF method

    We used the ‘‘SNFtool”in R software to achieve SNF,a useful and popular computational method for data integration in the field of disease subtype identification [19].It can deal with noise in different data types and make full use of common and complementary information across data types by integrating data in a non-linear way.We introduced SNF into multi-dimensional drug informatics data integration in this study for the first time,in the following three steps.1) DSNs are built for each data type.2) Multiple DSNs are integrated with SNF,and each of these DSNs is updated iteratively with information from other networks,making them more similar to each other than before.There are three main parameters in SNF:hyperparameter (η),number of neighbors (K),and number of iterations (T).The integration is robust to these parameters as described previously [19].Here we set η=0.5,K=20,and T=20,as recommended by Wang and colleagues [19].3) A final iDSN is obtained from SNF process convergence.

    Cluster validity index

    We used two cluster validity indexes to determine the number of clusters.

    Dunn index

    The Dunn index is the ratio of the smallest distance between observations not in the same cluster to the largest intracluster distance.The Dunn index has a value between zero and infinity and should be maximized.The Dunn index is calculated as follows:

    where C is the collection of all clusters,diam(Cm)is the largest intra-cluster distance in ClusterCm,whereasdist(Ck,Cl)is the distance between the nearest pair of samples in ClusterCkand ClusterCl.

    Silhouette index

    The silhouette value is a measure of how well each object lies within its cluster.The silhouette value ranges from-1 to 1,and should be maximized.It is calculated as follows:

    whereaiis the average distance between sampleiand all other data points within the same cluster,biis the lowest average distance of sampleito all points in any other clusters.

    Evaluation measurements

    Several evaluation measurements were used in the study,as introduced below.

    Normalized mutual information

    In probability theory and information theory,the mutual information (MI) of two random variables is a measure of their mutual dependence.In this study,Xis the ATC label vector of all the 593 drugs in iDSN,andYis the predicted label vector obtained by clustering these drugs.MI was calculated as follows:

    wherep(x,y)is the joint probability function ofXandY,whereasp(x)andp(y)are the marginal probability distribution functions ofXandY,respectively.Normalized MI(NMI) is calculated as follows:

    Connectivity

    The connectivity indicates the degree of connectedness of the clusters.DenoteNas the number of observations and denoteCas the collection of all clusters.Lrepresents the number of nearest neighbors.Definenni(j)as thejthnearest neighbor of observationi.Letxi,nni(j)be zero ifiandjare in the same cluster,and 1/jotherwise.The connectivity is defined as:

    The connectivity has a value between zero and infinity and should be minimized.

    Rogers-Tanimoto index

    The Rogers-Tanimoto similarity rely on a 2 × 2 contingency table,consisting of the following four cells:n11,n10,n01,andn00.n11is the number of observation pairs,where the two observations belong to the same cluster according to both partitionP1andP2.n10is the number of observation pairs,where the two observations belong to the same cluster according to partitionP1but not toP2.n01is the number of observation pairs,where the two observations belong to the same cluster according to partitionP2but not toP1.n00is the number of observation pairs,where the two observations do not belong to the same cluster according to both partitionP1andP2.The Rogers-Tanimoto similarity is defined as:

    Other integrative methods for comparison

    We compared the network fusion performance of PIMD with three previous integrative methods:1) the maximum method[21],2) the weighted average method [20,22,24-32],and 3)the probability disjunction [23].For the maximum method,we took the maximums of multiple drug similarity matrices.For the weighted average method,we averaged multiple drug similarity matrices by traversing weight.The weight of each drug similarity network is from 0 to 1 with step 0.1.For the probability disjunction,the formula is:

    whereSd,d′is the integrative similarity measurement between drugdandd′.

    Data type contribution

    For each edge in the iDSN,we used similarity scores from each single network to describe which data type was the primary contributor.First,we ranked three similarity scores of the edge in each single network asSi≥Sj≥Sk,wherei,j,andkrefer to the three types of data.IfSiwas ≥10% higher thanSj,the edge was attributed to theidata type.IfSiwas <10%higher thanSjbutSjwas ≥10%higher thanSk,the edge was attributed to both theiandjdata types.IfSiis <10%higher thanSj,butSjis <10%higher thanSk,the edge was attributed to all three data types.

    Statistical analysis

    We performed five types of enrichment analyses to annotate the drug clusters for drug precision repositioning.For drugbased and target-based enrichment analyses,we calculated enrichment score (ES) as follows:

    wherekis the number of drugs with a particular label (e.g.,ATC code) in the cluster of interest,mis the number of drugs with the label in the overall dataset,nis the total number of drugs in the cluster of interest,andNis the total number of drugs in the overall dataset.Then we used the hypergeometric distribution to calculate thePvalue as follows:

    Drug-based enrichment analysis

    For each cluster,we performed drug class and absorption,distribution,metabolism,excretion,and toxicity(ADMET)property enrichment analyses based on the DrugBank database.We computed ESs andPvalues based on drug ATC code,superclass label,and ADMET properties.

    Target-based enrichment analysis

    For target proteins of each cluster,we performed target class,Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway,and Gene Ontology (GO) enrichment analyses based on the International Union of Basic and Clinical Pharmacology/British Pharmacological Society (IUPHAR/BPS),KEGG,and GO databases,respectively[47-49].KEGG and GO analyses were conducted with the ClusterProfiler tool in R software[50].KEGG and GO enrichment results were selected by aPvalue threshold of 0.05.

    Drug property analysis

    Physicochemical features were extracted from a drug property list in RepurposeDB[51].Chemical descriptors were computed based on information from the Pybel,JOELib2,and Chemminer chemoinformatics libraries [45,52,53].We usedt-tests to evaluate the significance of deviation of mean values for each cluster from that of all drugs in DrugBank.In each cluster,for propertyiin property setI,the deviation from the mean value of all drugs in DrugBank is calculated as follows:

    whereA(i)is the mean property value for all drugs in Drug-Bank,andC(i)is the mean property value fro drugs in the cluster.

    Chemogenomic enrichment analysis

    We use the CGEA online tool to analyze drugs in each cluster(http://server.dudleylab.org/index).CGEA maps the drug list to various annotation resources,including drug-induced transcriptional modules,enzymes,and fragments.These potentially relevant features identified provide a rich chemogenomic context for drugs of interest.

    Chemical ontology enrichment analysis

    BiNChE is a web tool for chemical enrichment analysis based on the Chemical Entities of Biological Interest(ChEBI)Ontology [54].We performed enrichment analysis by mapping the drugs in each cluster to both ‘structure’ and ‘role’ subsets of ChEBI Ontology.

    Results

    An overview of PIMD framework

    PIMD is a systematic framework that predicts new therapeutic properties of known drugs,calculates the contributions of various types of data,and annotates various aspects of drug grouping,by integrating multiple drug properties.In this study,we hypothesize that if drugdAand drugdBare similar with respect to a particular parameter,then the therapeutic propertytof drugdAmay also be shared by drugdB.The functions of a drug can be characterized by multiple drug informatics.Therefore,in PIMD,we investigated three representative drug properties:1) chemical properties,based on chemical structure data derived from the PubChem Compound database[7];2)pharmacological properties,based on protein target sequence data derived from the DrugBank and UniProt databases [5,44];and 3) clinical and phenotypic properties,based on side effect data derived from the SIDER database [8].

    The PIMD consists of three parts (Figure 1).In the first part,we constructed three DSNs based on side effects,chemical structure,and drug targets,denoted as DSN-S,DSN-C,and DSN-T,respectively.Then,we applied the SNF method to fuse the DSN-S,DSN-C,and DSN-T into the iDSN in a manner that enables the relative contribution of each data type to be calculated.In the second part,we assessed drug repositioning using two approaches:finding drug pairs with a high similarity score,and identifying unexpected drugs in each cluster.In the iDSN,drug pairs are ranked according to their similarity scores,with highly ranked drug pairs being the most likely to achieve drug repositioning,wherein one drug in a pair may be repurposed for the therapeutic properties of the other drug.We used spectral clustering to distinguish clusters within the iDSN [55],wherein drugs within the same cluster have similar therapeutic properties, as reflected by first-level ATC codes. Within each cluster, drugs with unexpected ATC labels are flagged as having the potential for a repositioning of therapeutic properties. In the third part, aimed at drug precision repositioning, we carried out the following series of enrichment analyses to label and annotate each cluster: 1) drug-based enrichment analysis; 2) target-based enrichment analysis; 3)drug property analysis; 4) CGEA; and 5) chemical ontology enrichment analysis. These analyses can provide guidance to researchers conducting drug repositioning studies from various perspectives. They may also help to elucidate differences in the mode of action of different drug clusters and reveal potential associations between drugs within a cluster.

    Figure 1 Construction,drug repositioning,and annotation of iDSN using PIMD

    Global analysis of iDSN

    We used PubChem CIDs as the sole drug identifiers and extracted chemical structures, target sequences, and side effect sources for 593 drugs from the PubChem Compound, Drug-Bank, UniProt, and SIDER databases. After constructing the DSN-C, DSN-T, and DSN-S, and combining these three networks into the iDSN with the SNF method, we applied spectral clustering to divide the iDSN into 32 subnetworks (labeled Cluster 1-32) based on two validity indexes (Figure 2,Figure S1) [55,56].

    Figure 2 The iDSN and contributions of various data types

    Figure 3 Performance of the iDSN

    In the iDSN, each node represents a drug, and the edges connecting the nodes are thickness-weighted according to similarity of the connected drugs.Node color represents the cluster which the drug belongs to.Note that high-similarity connections are found predominantly within clusters.Edge color indicates the data type that is the main contributor to the similarity between the drugs,corresponding to the color scheme in the pie chart.Note that the iDSN model as a whole is supported by all analyzed types of data;in particular,drug pair similarities in the iDSN are supported by two or more types of data,and among the three single drug properties,the highest relative contribution comes from the drug target-based data.

    We divided the iDSN edges into two categories:withincluster edges and between-cluster edges.Data contribution analysis for edges within the cluster(Figure S2)shows that side effect-based data,chemical structure-based data,and drug target-based data account for 11.0%,22.2%,and 43.9%,respectively.Although the contributions of drug target-based data and side effect-based data are greater than that of chemical structure-based data for all edges in the iDSN (Figure 2),side effect-based data contribute the least to the edges inside the cluster and excessive contribution is from drug targetbased data (Figure S2).These results indicate that drug target-based data play a more important role in drug clustering.Data contribution details for each cluster are reported in Table S1.For instance,the contribution of chemical structure-based data in Cluster 6 accounts for 94.4%.We found that Cluster 6 drugs are all peptide drugs,which have chemical structures significantly different from those of other drugs.Cluster 6 drugs also exhibit a significant deviation for physicochemical features (Figure S3A and B).The deviation is consistent with the nature of these peptide drugs.

    Comparing the performances of the iDSN and each singleproperty network,we found that the iDSN has a clearer cluster structure than any of the three single-property DSNs(Figure 3A).To verify that,we calculated Dunn index,a metric for evaluating the quality of clustering results.The Dunn index for the iDSN is higher than those for the single-property DSNs(Figure 3B).

    Furthermore,compared to the single-property DSNs,the iDSN has a larger overlap between the drug cluster and drug ATC labels.NMI,which reflects consistency across the original ATC and cluster labels,was examined,and a higher NMI for the iDSN was obtained than those for the DSN-C,DSN-T,or DSN-S (Figure 3C).In addition,we compared the AOR of each single-property DSN and the iDSN by first ranking similarity scores of all drug pairs in the iDSN or single-property DSNs and assigning overlapping consecutive drug pair bins with 3000 pairs per bin,such that bin 1 contains the top 3000 most similar drug-drug pairs,bin 2 contains pairs ranked 100th to 3100th in similarity,and so on (first 320 bins are shown in Figure 3D).Subsequent calculation of AORs for each bin shows clearly that the AORs for the iDSN are higher than those of the single-property DSNs (Figure 3D).These results suggest that PIMD makes full use of common and complementary information about drug properties and that the iDSN performs better than any single-property DSNs.

    The network fusion of PIMD goes beyond a simple integration representing a maximum or an average of drug similarity measurements.It can capture potential links between drugs.On the one hand,if the similarity score in a single-property DSN is high,but the similarity scores in other singleproperty DSNs are low,PIMD does not dilute the original information.On the other hand,if the similarity scores for a drug pair in each of the three single-property DSNs are unremarkable,PIMD can still capture potential similarities.For example,the first-level ATC codes for fenoprofen and sulfasalazine are M and A,respectively.Fenoprofen (CID:000003342) is used for symptomatic relief of rheumatoid arthritis,osteoarthritis,and mild to moderate pain,whereas sulfasalazine (CID:005359476) is used to treat inflammatory bowel disease.The similarity scores for fenoprofen and sulfasalazine based on the three individual properties are relatively low (rank of 16,951 in DSN-C,13,466 in DSN-T,and 146,145 in DSN-S),while the iDSN similarity score ranked much higher (rank of 2729) and both drugs were placed in Cluster 8.Although crossover in an application would not have been predicted by any of the single-property DSNs,it has been reported that sulfasalazine can be used to treat rheumatoid arthritis [57,58].

    PIMD shows better performance than previous integrated methods

    To better evaluate the network fusion performance of PIMD,we compared our results with three previous integrative methods:1) the maximum method [21],2) the weighted average method [20,22,24-32],and 3) the probability disjunction [23].Here,we calculated two internal indices (connectivity and Dunn index) and two external comparison indices (NMI and Rogers-Tanimoto index).The internal indices are used to measure the goodness of a clustering structure without external information [59].The external indices are a measure of agreement between two partitions where the first partition is thea prioriknown clustering structure,and the second partition results from the clustering procedure [60].The connectivity indicates the degree of connectedness of the clusters.The Rogers-Tanimoto index reflects similarity between the original ATC and cluster labels.Among the four indices,only connectivity should be minimized.We found that PIMD performs the best compared with other methods (Figure 3E).The superior network fusion performance of PIMD results from the application of network-based approach,which is non-linear and utilizes topology information of the network.

    Drug repositioning from two aspects

    There are two approaches for drug repositioning:the first is finding drug pairs with a high similarity score;the second is finding unexpected drugs in each cluster.Drug pairs with high similarity in the iDSN could provide us with clues for drug repositioning.The AORs of the top 10 (Table 1) and the top 100 (top 1000 drug pairs are listed in Table S2) drug pairs reached 80%and 94%,respectively.Drug pairs with high similarity scores but different ATC codes may have the potential for repositioning.For example,triptorelin (CID:025074470)and nafarelin acetate (CID:025077649) are both gonadotropin-releasing hormone receptor agonists despite having different first-level ATC codes.Interestingly,this drug pair has a high iDSN similarity score and is ranked second in the drug association list.Within the top 20 drug pairs,7 drugs have been repurposed successfully according to RepurposeDB[51] and Repurposed Drug Database (http://drugrepurposingportal.com/),which record all repurposed drugs thus far.We further checked whether the repositioning in the databases is relevant to our prediction and found that 6 out of these 7 drugs are repositioned for the same purpose as we predicted.

    Unexpected drugs within an iDSN drug cluster can also signify the potential for repositioning based on distinct therapeutic properties.That is,a drug with an unusual ATC code within a certain cluster may be repositioned for alternative therapeutic properties.The unexpected drugs in each cluster are indicated in Table S3,and some cases are discussed in the section of drug repositioning case using PIMD.

    Table 1 Top 10 drug pairs in the iDSN

    We compared the two drug repositioning approaches together.In the top 100 drug pairs with a high similarity score,there are 6 drug pairs with totally different first-level ATC codes,among which 3 pairs can be discovered using the secondapproach.In the top 1000 drug pairs,there are 234 drug pairs with totally different first-level ATC codes,among which 75 pairs can be discovered using the second approach.These data indicate a certain overlap between the results predicted using these two approaches.The unexpected drugs in each cluster can achieve drug repositioning between different ATC codes.Compared to the second approach,top-ranking drug pairs with the same ATC codes are more inclined to achieve drug repositioning.

    To explore the global repositioning association among ATC codes,we analyzed the top 5%of drug pairs in the iDSN.After dividing these drugs into 14 groups according to their ATC codes,we averaged the similarity scores between the different drug groups and used these averages as indices of the repositioning potential within each ATC code.As shown inFigure 4(where edge thickness represents repositioning potential between two ATCs),there is a particularly high repositioning potential between dermatological drugs and respiratory system drugs.Among 78 dermatological drugs (D) and respiratory system drugs (R),23 drugs have been repurposed successfully between D and R (Table S4).

    Drug cluster annotation in the iDSN

    To label and annotate clusters for drug precision repositioning,we performed iDSN cluster analysis consisting of five statistical analysis methods from drug and target perspectives,thus providing multiple views to verify and select drug clusters for researchers in different fields.These analyses reveal potential links between drugs in the same clusters and differences in the modes of action between drugs in different clusters.

    ATC,superclass,and ADMET property enrichment analyses of drugs

    We set to explore to what extent drugs within a cluster share a common ATC code.For each cluster,we computed an ATC code enrichment score and an accompanyingPvalue.We found that 30 out of 32 clusters were significantly enriched for at least one ATC code (P<0.05;Figure 5A).Furthermore,25 out of 32 clusters were significantly enriched for at least one superclass label (P<0.05;Figure 5B).The superclass label of the drug is extracted from the DrugBank database and focuses more on the chemical attributes of drugs.Combining the two enrichment results may lead to new discoveries.For example,Cluster 22 was found to be enriched in the lignan/norlignan superclass as well as in the antineoplastic and immunomodulating agent ATC code.Recent studies have shown that lignans/norlignans play an important role in anticancer therapies [61].In addition,we extracted 18 ADMET property terms from the DrugBank database and examined whether drugs within the same cluster tend to have the same ADMET properties (Figure S3C).For example,6 clusters are enriched in the nervous system ATC code,while four of them (Cluster 11,Cluster 26,Cluster 28,and Cluster 30) are also enriched in the blood brain barrier(+).This observation is in line with the knowledge that nervous system drugs usually need to penetrate the blood brain barrier.

    Figure 4 Global repositioning association among ATC codes

    KEGG and GO enrichment analyses of targets

    To explore whether drugs in the same cluster tend to target similar proteins and whether particular drug classes are associated with particular target classes,we performed KEGG pathway and GO enrichment analyses for targets of each drug cluster [48-50].The most enriched pathways,biological processes,cellular components,and molecular functions(Figure 5C and D,Figure S3D and E;result matrices are listed in Table S5) differ substantially among the different clusters,highlighting differences in mode of action of drugs.Nonetheless,some drug clusters were found to contain common pathways and biological processes.We applied key biological themes for major clusters.For example,Cluster 27,which is enriched for antineoplastic and immunomodulating agents,was most enriched for the KEGG pathways of Rap1 signaling pathway,Ras signaling pathway,and central carbon metabolism in cancer.The most enriched GO biological process terms for Cluster 27 include protein autophosphorylation,positive regulation of MAPK cascade,and phosphatidylinositolmediated signaling.These pathways and biological processes are indeed closely related to the mechanisms of antineoplastic drugs.

    Another example is Cluster 16.Drugs in Cluster 16 are enriched for the ‘respiratory system’ ATC code.Chemical ontology enrichment analysis shows that drugs in this cluster are enriched for methylxanthine.The methylxanthines in Cluster 16,such as caffeine and theophylline,are used in therapy for respiratory diseases.However,the KEGG pathway enrichment analysis shows that drug targets in this cluster are also enriched for some cancer and immunology related pathways.These data suggest the association between methylxanthine and cancer,which is supported by recent studies [62,63].

    Furthermore,we conducted target enrichment analysis for each cluster based on class information.The class information is collected from the IUPHAR/BPS database [47],which indicates the type of drug targets,such as G-protein-coupled receptors (GPCR),catalytic receptor,and enzyme.We found that 24 of the 32 iDSN clusters were significantly enriched(P<0.05;Figure S3F).For example,Cluster 32 is enriched in the voltage-gated ion channel (VGIC) target class as well as in the ‘cardiovascular system’ ATC code.This observation is in line with the fact that most of Cluster 32 drugs are calcium channel blockers and used as antihypertensive drugs,whose targets mainly include calcium voltage-gated channel alpha1(CACNA1) subunits.

    Physicochemical feature and chemical descriptor analyses of drug properties

    Given pharmacological profiling of small molecules may also affect drug repositioning,we analyzed the characteristic physicochemical features and chemical descriptors for each cluster.In total,we extracted 14 physicochemical features from the RepurposeDB drug property list[51]and 62 chemical(i.e.,atomic,compositional,and geometric) descriptors.These properties were quantitated using the Pybel,JOELib2,and Chemminer chemoinformatics tools[45,52,53].For each property,we calculated mean values of these drug properties for all drugs in DrugBank and in each cluster.We found that the mean value of drug properties for drugs in each cluster was deviated from that for all drugs in DrugBank (Table S6).The degrees of statistically significant deviation (P<0.05;see Method) are shown in Figure S3A and G.Considering the bias resulting from the incompleteness of drug set,we also compared the mean value of drugs in each cluster with those of all the 593 drugs in the iDSN (Figure S3B and H).

    Figure 5 Drug-based and target-based enrichment analyses

    Chemogenomic enrichment analysis

    CGEA,similar to gene set enrichment analysis,is a method that compares drugs with a range of biological and chemical annotations[64](http://server.dudleylab.org/index).By identifying chemogenomic characteristics shared by sets of drugs such as enzymes,transporters,and structural fragments,we can obtain abundant chemogenomic context for a biological state of the drug set.We used CGEA to analyze drugs in each cluster for biological and chemical annotations.The results are listed in Table S7.The example of CGEA is described in detail in the section of ‘Case study:drug repositioning using 4 clusters of iDSN’.

    Enrichment analysis of chemical ontology

    Chemical ontology enrichment analysis is based on the ChEBI Ontology,which is a dictionary of chemical compounds with biological roles [65].Mapping drugs to the ChEBI database can improve our understanding of biochemical nature of drugs.We performed a chemical ontology enrichment analysis of the drugs in each cluster using the BiNChE tool [54].The resultant ChEBI-based enriched structure terms and role terms for each cluster are provided in Table S8.The example of chemical ontology enrichment analysis is described in detail in the section of ‘Case study:drug repositioning using 4 clusters of iDSN’.

    These analyses validate the rationality of the cluster division.Furthermore,by performing the 5 types of analyses described above,we have a comprehensive understanding of various properties of drugs in the cluster to conduct drug repositioning better.For example,Cluster 15 and Cluster 27 are both enriched for antineoplastic and immunomodulating agents,but their enriched KEGG pathways and GO terms(Figure 5C and D,Figure S3D and E) differed from each other.Chemical ontology enrichment analysis shows that Cluster 15 drugs have a role term of antimetabolite,while Cluster 27 drugs have a role term of protein kinase inhibitor.Most of Cluster 15 drugs are purine or pyrimidine analogs,whereas most of Cluster 27 drugs are tyrosine kinase inhibitors.These analyses highlight the differences in mode of action of the drugs.In Cluster 15,ribavirin (ATC code for ‘a(chǎn)ntiinfectives for systemic use’),which is a guanosine analog used for anti-virus,shows its therapeutic potential for cancers [66],and its mode of action is similar to other Cluster 15 drugs.

    Additionally,these analyses provide us a chance to look at the same drug cluster in combination of different perspectives,which may bring us some new discoveries.For example,Cluster 5 drugs are enriched for the ‘nervous system’ ATC code.Chemical ontology enrichment analysis shows that Cluster 5 drugs have a role term of serotonergic drug,a type of nervous system-related drugs.KEGG pathway enrichment and GO enrichment analyses also reveal terms related to the nervous system.However,there are also terms related to cardiovascular system according to cellular components (CC) and molecular function (MF) of GO enrichment analysis results(Figure S3D and E).In ADMET property enrichment analysis,Cluster 5 drugs are enriched in human ether-a-go-gorelated gene (hERG) inhibition (predictor I,strong inhibitor)and hERG inhibition (predictor II,inhibitor).hERG inhibition is related to QT prolongation.Moreover,CGEA results indicate that Cluster 5 drugs are enriched for the side effect electrocardiogram QT prolonged.These aforementioned analyses suggest that some serotonergic drugs(especially in Cluster 5) may cause QT prolongation.Indeed,some recent studies revealed the association between serotonergic drugs and QT prolongation [67,68].

    Case study:drug repositioning using 4 clusters of iDSN

    To illustrate the drug repositioning performance of PIMD,we examined some clusters in this section.

    The Cluster 28 drugs (Figure 2) are enriched for the ‘nervous system’ ATC code,the ‘benzenoids’ superclass label,and the ‘G-protein coupled receptors’ target class.KEGG pathway enrichment and GO enrichment analysis results also include terms related to the nervous system (Figure 6A and B,Figure S4A and B;Table S9).Our CGEA analysis shows that Cluster 28 drugs are enriched for the molecular fragment CCCN(C)C and drug-induced transcriptional module PC3-3(Figure 6C).The module contains expression profiles for 38 genes in response to 25 drugs,6 drugs out of which are included in Cluster 28 [69].The most prominent drug mode of action among Cluster 28 drugs is ‘Antihistamines for systemic use’,which can affect the central nervous system [70].Chemical ontology enrichment analysis illustrates that Cluster 28 drugs are enriched with neurotransmitters,neurotransmitter derivatives,and central nervous system drugs (Figure 6D-F;drug property analysis results in Figure 6F).We found that cyclobenzaprine,a drug used to treat skeletal muscle spasms and fibromyalgia,is an unexpected drug in Cluster 28.Cyclobenzaprine,as well as another drug named the antidepression drug amitriptyline in this cluster,exhibits antagonistic effects on the 5-hydroxytryptamine receptor 2A and possesses a tricyclic structure.Cyclobenzaprine is being studied for the post-traumatic stress disorder treatment according to ClinicalTrials.gov (https://clinicaltrials.gov/),a database of clinical studies conducted around the world.Therefore,the novel therapeutic property for cyclobenzaprine might be N(nervous system).These results indicate that Cluster 28 drugs are related to nervous system and show the good compatibility.Unexpected drugs in Cluster 28 also have the potential to be repositioned to ‘nervous system’ ATC code.

    Similarly,pentoxifylline,an unexpected drug in Cluster 16(drugs in this cluster are enriched for the ‘respiratory system’ATC code),is an interesting case.Pentoxifylline carries an ATC code for ‘cardiovascular system’ and is used in therapy for intermittent claudication [71].Notably,it has a high similarity score with theophylline (ATC code for ‘respiratory system’),another drug in Cluster 16 used to treat respiratory diseases such as asthma.They are both members of the xanthine family,so they have a relatively high similarity score based on chemical structure (rank of 2898 in DSN-C).They also have a relatively high similarity score based on drug targets (rank of 7277 in DSN-T) because of the presence of 5 common targets.Pentoxifylline can increase red blood cell deformability and decrease blood viscosity[72].But its similarity with theophylline has encouraged researchers to explore the potential of pentoxifylline to treat asthma [73].

    Clusters simultaneously enriched in two therapeutic properties can also provide a unique perspective for drug repurposing.Cluster 3 drugs (Figure 2) are enriched for both the‘cardiovascular system’ and ‘genitourinary system and reproductive hormones’ ATC codes (Figure S5A-G;Table S9).In this cluster,iloperidone,a drug used for schizophrenia,may be repositioned for applications related to the ‘cardiovascular system’ and the ‘genitourinary system and reproductive hormones’.Indeed,a previous study showed that repeat administration of iloperidone moderated hypotension [74].Interestingly,we noted that sildenafil,a drug that was successfully repositioned from ‘cardiovascular system’ treatment to‘genitourinary system and reproductive hormones’ treatment,is also present in Cluster 3,suggesting that drugs in this cluster may have the potential for both therapeutic properties.Likewise,prazosin and terazosin are an interesting pair.Prazosin(ATC code for ‘cardiovascular system’) is used to treat hypertension and also for urinary hesitancy associated with prostatic hyperplasia [75].Terazosin (ATC code for ‘genitourinary system and reproductive hormones’) is used to treat enlarged prostate symptoms,which can also moderate hypertension[76,77].Hence,the primary indications of these two drugs can be treated as secondary indications of each other.In terms of target proteins,both drugs have antagonistic effects on alpha-1A,-1B,and -1D adrenergic receptors.In terms of chemical structures,their most common substructure is the molecular fragment O=CN1CCNCC1.Their common adverse reactions include dizziness,headache,drowsiness,lack of energy,and weakness.

    Figure 6 Community analysis of Cluster 28

    Another interesting case is cyproterone acetate (ATC code for ‘genito-urinary system and sex hormones’) in Cluster 19.Cyproterone acetate can not only treat androgen-dependent conditions like excessive hair growth and acne but also treat prostate cancer[78].Therefore,the novel therapeutic property for cyproterone acetate could be L (antineoplastic and immunomodulating agents).More interestingly,cyproterone acetate was originally developed as a progestin [79],but it was first marketed as an antiandrogen [80].Cluster 19 drugs are enriched for both the ‘genito-urinary system and sex hormones’ and ‘a(chǎn)ntineoplastic and immunomodulating agents’ATC codes.In this cluster,cyproterone acetate has high similarity scores with bicalutamide (ATC code for ‘a(chǎn)ntineoplastic and immunomodulating agents’) and medroxyprogesterone acetate(ATC code for‘a(chǎn)ntineoplastic and immunomodulating agents’ and ‘genito-urinary system and sex hormones’).For the drug pair of cyproterone acetate and medroxyprogesterone acetate,their similarity score based on the chemical structure is relatively high(rank of 101 in DSN-C),but the similarity score based on drug targets is relatively low(rank of 35,974 in DSNT).On the contrary,cyproterone acetate and bicalutamide have relatively high similarity score based on drug targets(rank of 178 in DSN-T) but relatively low similarity score based on chemical structures(rank of 81,257 in DSN-C).This is because cyproterone acetate is more structurally similar to medroxyprogesterone acetate (a type of progestin) but has antiandrogenic effects on androgen receptor like bicalutamide(a kind of antiandrogen).This explains the uniqueness of cyproterone acetate.These results also show that PIMD takes advantage of the universality,individuality,and complementarity of different drug properties.

    Discussion

    Here,we proposed PIMD,an integrative,systematic,and extensible framework for discovering novel therapeutic properties of drugs from heterogeneous data sources.PIMD characterizes the iDSN by integrating chemical structure data,target protein sequence data,and side effect data;and performing spectral clustering of the iDSN identified 32 drug clusters.Additionally,PIMD facilitates drug repositioning from two aspects:drug pairs with high iDSN similarity score and unexpected drugs in each cluster.Finally,via a series of enrichment analyses,PIMD annotates and evaluates all clusters from chemical,pharmacological,and genomic views.Thus,PIMD screens suitable clusters for drug precision repositioning.

    By integrating multi-dimensional drug informatics data,PIMD can capture potential similarities between drugs with sensitivity.The iDSN is superior to single-property DSNs in multiple evaluation measurements,including AOR,NMI,and cluster compactness.In this study,we primarily used ATC label as a golden standard to annotate drugs and evaluate the performance of PIMD and other compared methods.Interestingly,if the superclass label is used to calculate the evaluation measurements instead of ATC code,the SOR and NMI values of chemical structure-based DSN are higher than those of other DSNs,even the iDSN(Figure S6A and B).This suggests that chemical property has a large positive effect on superclass,whereas clinical and pharmacological properties could affect drug superclass label negatively.Indeed,ATC label and superclass label are two different drug catalogs.ATC label can comprehensively reflect therapeutic properties of drugs,thus commonly used to classify and label drugs,whereas superclass label primarily represents the chemical properties of the drugs.This result can improve our understanding of the difference between ATC and superclass labels of the drugs.

    PIMD offers new insights into the universality,individuality,and complementarity of three drug properties,including chemical,pharmacological,and clinical properties.Our calculations of the relative contribution of each data type indicate that the iDSN is driven by all data types.PIMD makes full use of the information on each property.Examination of NMI among single-property,dual-property,and threeintegrated-property DSNs (Figure S7A-D) shows that networks based on single property alone overlap marginally but are complementary to each other instead.Furthermore,PIMD provides the drug property contribution to each cluster,improving our understanding of the cluster characteristics(Table S1).

    Rather than simple fusion,PIMD provides a systematic way to evaluate the drug cluster and drug repositioning.Pharmacochemistry and pharmacogenomics researchers can use PIMD to screen drug clusters based on their own requirements(Figure 5A-D,Figure S3A-F;Tables S7 and S8).

    As a highly extensible framework,PIMD can fuse various properties beyond the three properties examined here.For example,combining drug expression profiles would allow us to elucidate similarities between drugs at the transcriptional level.Drugs have multi-dimensional properties,and drug effects on disease processes are an interdisciplinary issue.The complementarity among multiple properties allows us to assess drug similarities more accurately and thus to provide a more comprehensive and clearer direction for drug repositioning.

    There are two reasons why we chose these three data sources as an example.Firstly,these three data sources are representative of a variety of drug informatics data.Side effect,chemical structure,and molecular target data represent the clinical,chemical,and pharmacological properties of drugs,respectively.The three properties comprehensively summarize the drug characteristics.Secondly,these three data sources have sufficient data to extract.We can collect plenty of these data from the public databases.Moreover,despite insufficient data available,we also compared iDSN with DSNs based on three other data types:1)drug 3D structure,2)drug expression profiles,and 3) PPI network (similarity measurements can be found in File S1).The results show that iDSN performs best compared with other DSNs (Figure S8).

    The current study is limited by the set of drugs available when using the intersection of drug sets with multidimensional properties.If a certain property of a drug is not available,it would be excluded from the construction of iDSN.Therefore,with the accumulation of drug informatics data in the future,we expect that the scale of iDSN would be expanded,and the performance of PIMD would be further improved accordingly.

    In summary,PIMD provides a new perspective for drug repositioning through multi-property fusion and an analysis package.It facilitates to understand the integration of drug properties at a deeper level,and its high expansibility and modularity would allow users to explore drugs from a wider range of fields.

    Code availability

    Source code of PIMD is available at https://github.com/Sepstar/PIMD/.

    CRediT author statement

    Song He:Conceptualization,Methodology,Investigation,Software,Visualization,Writing -original draft,Writing -review& editing.Yuqi Wen:Conceptualization,Methodology,Investigation,Software,Visualization,Writing-original draft,Writing-review& editing.Xiaoxi Yang:Investigation,Visualization.Zhen Liu:Investigation,Visualization.Xinyu Song:Investigation,Visualization.Xin Huang:Investigation,Visualization.Xiaochen Bo:Conceptualization,Methodology,Supervision,Software,Writing -review & editing.All authors read and approved the final manuscript.

    Competing interests

    The authors have declared no competing interests.

    Acknowledgments

    This work was supported by the National Natural Science Foundation of China(Grant No.U1435222)and the Program of International Sci-Tech Cooperation,China (Grant No.2014DFB30020).

    Supplementary material

    Supplementary data to this article can be found online at https://doi.org/10.1016/j.gpb.2018.10.012.

    ORCID

    0000-0002-4136-6151 (Song He)

    0000-0002-2391-3297 (Yuqi Wen)

    0000-0002-2674-1902 (Xiaoxi Yang)

    0000-0002-3327-4118 (Zhen Liu)

    0000-0002-1905-3503 (Xinyu Song)

    0000-0002-6569-0095 (Xin Huang)

    0000-0003-1911-7922 (Xiaochen Bo)

    999久久久精品免费观看国产| 国产伦在线观看视频一区| 两个人视频免费观看高清| 国产1区2区3区精品| 一本一本综合久久| 精品久久久久久久毛片微露脸| 中文字幕高清在线视频| 香蕉国产在线看| 久久久水蜜桃国产精品网| 日韩成人在线观看一区二区三区| 精品少妇一区二区三区视频日本电影| 国产av一区在线观看免费| 国产精品精品国产色婷婷| 九九热线精品视视频播放| 午夜福利欧美成人| 一进一出好大好爽视频| 亚洲专区国产一区二区| 黑人巨大精品欧美一区二区mp4| 亚洲成人久久性| 亚洲一区二区三区不卡视频| 成人一区二区视频在线观看| 亚洲国产中文字幕在线视频| 国产麻豆成人av免费视频| 18美女黄网站色大片免费观看| 久久亚洲真实| 国产成人精品久久二区二区91| 国产黄色小视频在线观看| 丰满的人妻完整版| 身体一侧抽搐| 人人妻人人看人人澡| 欧美中文综合在线视频| 久久九九热精品免费| 欧美成人性av电影在线观看| 女生性感内裤真人,穿戴方法视频| 琪琪午夜伦伦电影理论片6080| 免费电影在线观看免费观看| 精品乱码久久久久久99久播| 99久久精品热视频| av在线播放免费不卡| 国产一级毛片七仙女欲春2| 午夜激情av网站| 国产精品野战在线观看| 欧美一级a爱片免费观看看 | 波多野结衣巨乳人妻| 一本综合久久免费| 熟女电影av网| 国产精品亚洲一级av第二区| 美女扒开内裤让男人捅视频| 嫁个100分男人电影在线观看| 九九热线精品视视频播放| 欧美一区二区国产精品久久精品 | 欧美乱妇无乱码| 国产欧美日韩一区二区精品| xxxwww97欧美| www.熟女人妻精品国产| 欧美大码av| tocl精华| 免费看美女性在线毛片视频| 搡老妇女老女人老熟妇| 中文字幕熟女人妻在线| videosex国产| 久久久久久大精品| 久久人妻福利社区极品人妻图片| а√天堂www在线а√下载| 免费搜索国产男女视频| 欧美成狂野欧美在线观看| 久久久精品国产亚洲av高清涩受| 十八禁人妻一区二区| 国产精品精品国产色婷婷| 国产私拍福利视频在线观看| av欧美777| 国产精品av视频在线免费观看| 丰满人妻一区二区三区视频av | 老司机靠b影院| 一进一出抽搐gif免费好疼| 亚洲精品在线美女| 老司机深夜福利视频在线观看| 法律面前人人平等表现在哪些方面| 首页视频小说图片口味搜索| 在线播放国产精品三级| 免费在线观看日本一区| 黑人欧美特级aaaaaa片| 日日夜夜操网爽| 精品不卡国产一区二区三区| 久久久水蜜桃国产精品网| 99精品久久久久人妻精品| 免费在线观看日本一区| 男女午夜视频在线观看| 夜夜爽天天搞| 看免费av毛片| 精品国产超薄肉色丝袜足j| 午夜福利18| 草草在线视频免费看| 老鸭窝网址在线观看| 亚洲国产精品久久男人天堂| 久久人妻av系列| 国产激情偷乱视频一区二区| 窝窝影院91人妻| 男人的好看免费观看在线视频 | 成人国语在线视频| 成年女人毛片免费观看观看9| 校园春色视频在线观看| 国产精品亚洲一级av第二区| 19禁男女啪啪无遮挡网站| 国产成人aa在线观看| 精品第一国产精品| 日韩精品中文字幕看吧| 小说图片视频综合网站| 亚洲欧美精品综合一区二区三区| 国产成+人综合+亚洲专区| 两个人看的免费小视频| 男女那种视频在线观看| 久久精品成人免费网站| 亚洲av成人精品一区久久| 国产一区二区三区视频了| 99国产极品粉嫩在线观看| 草草在线视频免费看| 91麻豆精品激情在线观看国产| 少妇人妻一区二区三区视频| 12—13女人毛片做爰片一| 蜜桃久久精品国产亚洲av| 男女视频在线观看网站免费 | 国产三级黄色录像| 一级a爱片免费观看的视频| 国产一区二区三区视频了| 一边摸一边做爽爽视频免费| 美女免费视频网站| 丰满的人妻完整版| 视频区欧美日本亚洲| 亚洲男人天堂网一区| 婷婷精品国产亚洲av在线| 日本一本二区三区精品| 最近最新免费中文字幕在线| 一a级毛片在线观看| 欧美黑人欧美精品刺激| 天堂√8在线中文| 极品教师在线免费播放| 亚洲av电影不卡..在线观看| 全区人妻精品视频| 日本精品一区二区三区蜜桃| 最新在线观看一区二区三区| 麻豆成人av在线观看| 国产精品,欧美在线| 亚洲av熟女| 美女 人体艺术 gogo| 国产精品久久电影中文字幕| 无人区码免费观看不卡| 麻豆成人午夜福利视频| 1024手机看黄色片| x7x7x7水蜜桃| 成人18禁在线播放| 亚洲成人免费电影在线观看| 淫秽高清视频在线观看| 日韩欧美在线二视频| 成人特级黄色片久久久久久久| 国内揄拍国产精品人妻在线| 亚洲狠狠婷婷综合久久图片| 精品久久久久久久末码| 欧美成人一区二区免费高清观看 | 久久久久亚洲av毛片大全| 精品人妻1区二区| 亚洲,欧美精品.| 又大又爽又粗| 女人爽到高潮嗷嗷叫在线视频| 村上凉子中文字幕在线| 久久久久久大精品| 后天国语完整版免费观看| 国产成人av教育| 国产精品亚洲av一区麻豆| 91大片在线观看| 一级毛片精品| 亚洲成av人片在线播放无| 国产精品国产高清国产av| 久久国产乱子伦精品免费另类| 大型av网站在线播放| 久久午夜综合久久蜜桃| 亚洲色图av天堂| 午夜福利18| 国内揄拍国产精品人妻在线| 亚洲精品在线观看二区| 国产免费男女视频| 伦理电影免费视频| 亚洲人成电影免费在线| 毛片女人毛片| 校园春色视频在线观看| 国产激情久久老熟女| 熟女少妇亚洲综合色aaa.| www.www免费av| 亚洲黑人精品在线| 男女之事视频高清在线观看| 亚洲av日韩精品久久久久久密| 国产一区在线观看成人免费| 曰老女人黄片| 国内精品久久久久精免费| 美女扒开内裤让男人捅视频| 色综合欧美亚洲国产小说| www国产在线视频色| 亚洲av成人一区二区三| 1024香蕉在线观看| 别揉我奶头~嗯~啊~动态视频| 亚洲人成网站高清观看| 国产精品亚洲美女久久久| 亚洲一区中文字幕在线| 特级一级黄色大片| 亚洲专区字幕在线| 欧美日韩黄片免| 五月伊人婷婷丁香| 国产精品一区二区三区四区久久| 欧美黑人巨大hd| 国产成人av激情在线播放| 日韩欧美 国产精品| 少妇被粗大的猛进出69影院| 亚洲人成77777在线视频| 99国产综合亚洲精品| 最近最新中文字幕大全免费视频| 免费观看人在逋| 99riav亚洲国产免费| 欧美日韩瑟瑟在线播放| 12—13女人毛片做爰片一| 日本精品一区二区三区蜜桃| videosex国产| 久久伊人香网站| 99国产精品一区二区三区| 亚洲国产中文字幕在线视频| 国产精品 欧美亚洲| 成年免费大片在线观看| 欧美一级毛片孕妇| 啦啦啦观看免费观看视频高清| 不卡av一区二区三区| 亚洲电影在线观看av| 舔av片在线| 精品久久久久久成人av| 国产精品一区二区精品视频观看| 日韩成人在线观看一区二区三区| 1024香蕉在线观看| 亚洲国产欧洲综合997久久,| 免费看美女性在线毛片视频| 亚洲欧美一区二区三区黑人| 真人一进一出gif抽搐免费| 国产精品亚洲av一区麻豆| 可以在线观看毛片的网站| 午夜福利18| 久久天躁狠狠躁夜夜2o2o| 精品国产乱子伦一区二区三区| 国产亚洲精品一区二区www| 久久精品91蜜桃| av天堂在线播放| 黄色女人牲交| av国产免费在线观看| 琪琪午夜伦伦电影理论片6080| 一夜夜www| 一进一出抽搐动态| 欧美激情久久久久久爽电影| 久久香蕉国产精品| 最近在线观看免费完整版| 国产又黄又爽又无遮挡在线| 精品第一国产精品| 亚洲成av人片在线播放无| 国产伦在线观看视频一区| 又黄又爽又免费观看的视频| 村上凉子中文字幕在线| 久久天堂一区二区三区四区| 国产高清有码在线观看视频 | 日本熟妇午夜| 国产午夜精品久久久久久| 欧美中文综合在线视频| 精品人妻1区二区| 18美女黄网站色大片免费观看| 丰满人妻熟妇乱又伦精品不卡| 老司机在亚洲福利影院| 999久久久国产精品视频| 男人舔女人的私密视频| 欧美日本亚洲视频在线播放| 男女下面进入的视频免费午夜| 国产精品一及| 欧美zozozo另类| 99热这里只有精品一区 | 亚洲中文字幕一区二区三区有码在线看 | 在线观看日韩欧美| 精品国产美女av久久久久小说| 国产精品九九99| 51午夜福利影视在线观看| 欧美日韩一级在线毛片| 欧美一级毛片孕妇| 中文资源天堂在线| 日韩成人在线观看一区二区三区| 欧美zozozo另类| 日本在线视频免费播放| 日本黄色视频三级网站网址| 久久婷婷人人爽人人干人人爱| 大型av网站在线播放| 亚洲av第一区精品v没综合| 午夜老司机福利片| 怎么达到女性高潮| 亚洲精品在线美女| 岛国视频午夜一区免费看| 一级片免费观看大全| 免费看日本二区| 婷婷精品国产亚洲av| 亚洲片人在线观看| 91成年电影在线观看| aaaaa片日本免费| 一本大道久久a久久精品| 国产精品av视频在线免费观看| 两性夫妻黄色片| 欧美国产日韩亚洲一区| 黄色成人免费大全| 国产高清激情床上av| 91成年电影在线观看| 午夜亚洲福利在线播放| 两性夫妻黄色片| 99热这里只有是精品50| 搞女人的毛片| 日韩 欧美 亚洲 中文字幕| 色精品久久人妻99蜜桃| 日本熟妇午夜| 黄色视频不卡| x7x7x7水蜜桃| 欧美日韩国产亚洲二区| 97人妻精品一区二区三区麻豆| 国产亚洲欧美在线一区二区| 中国美女看黄片| 国产成人av教育| 亚洲av片天天在线观看| 成人av一区二区三区在线看| 美女高潮喷水抽搐中文字幕| 欧美成人免费av一区二区三区| 国产午夜福利久久久久久| 久久精品综合一区二区三区| 欧美黄色淫秽网站| 免费观看精品视频网站| 国产熟女xx| 小说图片视频综合网站| 色噜噜av男人的天堂激情| 久久久水蜜桃国产精品网| 国产成人一区二区三区免费视频网站| 亚洲一区中文字幕在线| 男女做爰动态图高潮gif福利片| 黄色片一级片一级黄色片| 日日夜夜操网爽| 亚洲专区中文字幕在线| 国产一区二区三区在线臀色熟女| 亚洲精品久久国产高清桃花| 免费看日本二区| 人人妻人人看人人澡| 91字幕亚洲| 国产成年人精品一区二区| 欧美黑人巨大hd| 亚洲va日本ⅴa欧美va伊人久久| 90打野战视频偷拍视频| av片东京热男人的天堂| 国产精品免费视频内射| 午夜福利欧美成人| 欧美黑人欧美精品刺激| 亚洲av第一区精品v没综合| 亚洲一卡2卡3卡4卡5卡精品中文| 婷婷丁香在线五月| 亚洲欧美精品综合一区二区三区| 久久天堂一区二区三区四区| 伦理电影免费视频| 三级国产精品欧美在线观看 | 亚洲欧美精品综合一区二区三区| 国产精品久久久久久亚洲av鲁大| 男女床上黄色一级片免费看| 精品久久久久久久末码| 女警被强在线播放| 亚洲专区字幕在线| 久久久久久国产a免费观看| 久久伊人香网站| 91字幕亚洲| 制服丝袜大香蕉在线| 国内久久婷婷六月综合欲色啪| 手机成人av网站| 人人妻人人澡欧美一区二区| 国产精品av久久久久免费| 午夜老司机福利片| 久久午夜综合久久蜜桃| 又黄又爽又免费观看的视频| 舔av片在线| 久久久久精品国产欧美久久久| 天天躁夜夜躁狠狠躁躁| or卡值多少钱| 午夜免费激情av| 亚洲国产精品久久男人天堂| 成人国产综合亚洲| 国产亚洲精品第一综合不卡| 亚洲熟妇熟女久久| 无人区码免费观看不卡| 亚洲中文字幕日韩| 亚洲七黄色美女视频| 淫妇啪啪啪对白视频| 国产高清视频在线播放一区| 国产精品,欧美在线| 欧美成人性av电影在线观看| 首页视频小说图片口味搜索| 久久精品国产综合久久久| 国产高清视频在线播放一区| 久久人妻福利社区极品人妻图片| 婷婷丁香在线五月| 国产精品久久久av美女十八| 亚洲国产中文字幕在线视频| 中文亚洲av片在线观看爽| 听说在线观看完整版免费高清| 最近在线观看免费完整版| 色尼玛亚洲综合影院| 999久久久精品免费观看国产| 床上黄色一级片| 成人亚洲精品av一区二区| 成人特级黄色片久久久久久久| 88av欧美| 午夜久久久久精精品| 亚洲成人中文字幕在线播放| 香蕉国产在线看| 国产成+人综合+亚洲专区| 美女免费视频网站| 久久久久国内视频| 淫妇啪啪啪对白视频| 欧美三级亚洲精品| 亚洲性夜色夜夜综合| 成人三级做爰电影| 久热爱精品视频在线9| 我要搜黄色片| 超碰成人久久| 波多野结衣高清作品| 无人区码免费观看不卡| 国产成人精品久久二区二区免费| 97人妻精品一区二区三区麻豆| 国产亚洲精品久久久久5区| 非洲黑人性xxxx精品又粗又长| 国产精品国产高清国产av| 2021天堂中文幕一二区在线观| 精品国产乱码久久久久久男人| 最近最新中文字幕大全电影3| 99精品欧美一区二区三区四区| 美女 人体艺术 gogo| 欧美成人午夜精品| avwww免费| 久久性视频一级片| 99国产精品一区二区三区| 麻豆成人av在线观看| 精华霜和精华液先用哪个| 免费看美女性在线毛片视频| 日本三级黄在线观看| 999久久久精品免费观看国产| 欧美日韩福利视频一区二区| 国产三级中文精品| 天天添夜夜摸| 国产单亲对白刺激| 特级一级黄色大片| 久久中文看片网| 两性夫妻黄色片| 成人一区二区视频在线观看| 久久久久亚洲av毛片大全| 这个男人来自地球电影免费观看| 日日夜夜操网爽| 国产真人三级小视频在线观看| 欧洲精品卡2卡3卡4卡5卡区| 国产真人三级小视频在线观看| 国产久久久一区二区三区| 亚洲精品在线观看二区| 免费在线观看视频国产中文字幕亚洲| 丰满人妻一区二区三区视频av | 午夜福利在线观看吧| 欧美久久黑人一区二区| 久久久久久免费高清国产稀缺| 亚洲狠狠婷婷综合久久图片| 九色国产91popny在线| 美女免费视频网站| 国模一区二区三区四区视频 | 婷婷丁香在线五月| 国产主播在线观看一区二区| 日本一本二区三区精品| 国产黄a三级三级三级人| 亚洲av成人精品一区久久| 亚洲国产看品久久| 亚洲成a人片在线一区二区| 中文亚洲av片在线观看爽| 90打野战视频偷拍视频| 色综合欧美亚洲国产小说| 999精品在线视频| 国产三级黄色录像| 国内久久婷婷六月综合欲色啪| 女人高潮潮喷娇喘18禁视频| 99久久无色码亚洲精品果冻| 日韩欧美在线乱码| 中文字幕人妻丝袜一区二区| 久久精品国产综合久久久| 男女之事视频高清在线观看| 国产视频一区二区在线看| 成年人黄色毛片网站| 午夜免费成人在线视频| 亚洲欧美精品综合久久99| 丰满人妻一区二区三区视频av | 亚洲精品一区av在线观看| 久久久久国内视频| 国产精品国产高清国产av| 国产欧美日韩一区二区三| 女警被强在线播放| www日本在线高清视频| 色综合站精品国产| 成人18禁高潮啪啪吃奶动态图| 99re在线观看精品视频| 久9热在线精品视频| 欧美乱码精品一区二区三区| 在线观看午夜福利视频| 小说图片视频综合网站| 亚洲精华国产精华精| 国产乱人伦免费视频| 香蕉丝袜av| 亚洲av中文字字幕乱码综合| 美女午夜性视频免费| 亚洲 欧美一区二区三区| 国产精品亚洲美女久久久| 90打野战视频偷拍视频| 中文字幕久久专区| 国产蜜桃级精品一区二区三区| 亚洲一区高清亚洲精品| 欧美三级亚洲精品| 日韩高清综合在线| 午夜日韩欧美国产| 国产av不卡久久| 一卡2卡三卡四卡精品乱码亚洲| 人妻夜夜爽99麻豆av| 亚洲人成77777在线视频| 精品一区二区三区av网在线观看| 日韩 欧美 亚洲 中文字幕| 久久九九热精品免费| 亚洲精品中文字幕在线视频| 亚洲精品久久国产高清桃花| 好男人电影高清在线观看| av片东京热男人的天堂| 国产私拍福利视频在线观看| 国产成人av激情在线播放| 亚洲成人国产一区在线观看| 美女扒开内裤让男人捅视频| 成人三级黄色视频| 麻豆av在线久日| 国产精品电影一区二区三区| 丁香欧美五月| 黄频高清免费视频| 国产精品永久免费网站| 国产三级中文精品| 午夜精品在线福利| 悠悠久久av| 亚洲精品中文字幕一二三四区| 国产亚洲精品av在线| 亚洲精品粉嫩美女一区| 好看av亚洲va欧美ⅴa在| 一个人观看的视频www高清免费观看 | 好男人在线观看高清免费视频| 成人三级做爰电影| 男女午夜视频在线观看| 国产成人av教育| 精华霜和精华液先用哪个| 熟女少妇亚洲综合色aaa.| 日韩欧美一区二区三区在线观看| 色综合欧美亚洲国产小说| 大型黄色视频在线免费观看| 成年人黄色毛片网站| 亚洲欧美日韩高清在线视频| 亚洲熟妇熟女久久| 日本三级黄在线观看| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲人成电影免费在线| 妹子高潮喷水视频| 亚洲av美国av| 特大巨黑吊av在线直播| 亚洲av成人不卡在线观看播放网| 琪琪午夜伦伦电影理论片6080| 一级毛片女人18水好多| 每晚都被弄得嗷嗷叫到高潮| 久久性视频一级片| 小说图片视频综合网站| tocl精华| 一二三四在线观看免费中文在| 国产成人影院久久av| 国语自产精品视频在线第100页| 色精品久久人妻99蜜桃| 国产片内射在线| 一进一出好大好爽视频| 成人亚洲精品av一区二区| 久久精品国产清高在天天线| 欧洲精品卡2卡3卡4卡5卡区| 国产免费av片在线观看野外av| 国产久久久一区二区三区| 深夜精品福利| 精品一区二区三区av网在线观看| videosex国产| 久久精品亚洲精品国产色婷小说| 最近最新免费中文字幕在线| 久久久久久免费高清国产稀缺| 免费人成视频x8x8入口观看| 91九色精品人成在线观看| 在线永久观看黄色视频| 蜜桃久久精品国产亚洲av| 丁香欧美五月| 亚洲一区高清亚洲精品| 女人高潮潮喷娇喘18禁视频| 18美女黄网站色大片免费观看| www日本在线高清视频| 亚洲精品中文字幕一二三四区| 欧美日韩国产亚洲二区| 亚洲av日韩精品久久久久久密| 在线国产一区二区在线| 欧美性猛交╳xxx乱大交人| 国产精品野战在线观看| 亚洲乱码一区二区免费版| 一级作爱视频免费观看| 精品一区二区三区视频在线观看免费| 无人区码免费观看不卡| 色在线成人网| 成人18禁在线播放| 老司机午夜十八禁免费视频| 999久久久国产精品视频| av福利片在线|