
==== Front
Sci Data
Sci Data
Scientific Data
2052-4463
Nature Publishing Group UK London

3734
10.1038/s41597-024-03734-x
Data Descriptor
KORE-Map 1.0: Korean medicine Omics Resource Extension Map on transcriptome data of tonifying herbal medicine
http://orcid.org/0000-0002-5383-5285
Park Musun 1
http://orcid.org/0000-0003-2915-0742
Park Sang-Min 2
Lee Haeseung 3
Kim Aeyung 4
Kim No Soo 5
Kim Yu Ri 5
Yi Jin-Mu jmyi@kiom.re.kr

5
Cha Seongwon scha@kiom.re.kr

1
1 https://ror.org/005rpmt10 grid.418980.c 0000 0000 8749 5149 Korean Medicine (KM) Data Division, Korea Institute of Oriental Medicine, Daejeon, Republic of Korea
2 https://ror.org/0227as991 grid.254230.2 0000 0001 0722 6377 College of Pharmacy, Chungnam National University, Daejeon, Republic of Korea
3 https://ror.org/01an57a31 grid.262229.f 0000 0001 0719 8572 College of Pharmacy, Pusan National University, Busan, Republic of Korea
4 https://ror.org/005rpmt10 grid.418980.c 0000 0000 8749 5149 KM Application Center, Korea Institute of Oriental Medicine, Daegu, Republic of Korea
5 https://ror.org/005rpmt10 grid.418980.c 0000 0000 8749 5149 KM Convergence Research Division, Korea Institute of Oriental Medicine, Daejeon, Republic of Korea
6 9 2024
6 9 2024
2024
11 97420 3 2024
2 8 2024
© The Author(s) 2024
2024
https://creativecommons.org/licenses/by-nc-nd/4.0/ Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.
Traditional herbal medicine, rooted in a long history of use in East Asia, combines several herbs to create treatments showing high efficacy with minimal side effects, for specific diseases. Such combination therapies represent a potential reservoir of new drugs for treating multifactorial and incurable chronic diseases. However, the complexity of their mechanisms of action due to the combination of multiple compounds, has limited their research integration into modern pharmacological science. To address this challenge, we constructed drug-induced transcriptome data for herbal medicines through systematic experiments, analyzed with the aid of various omics databases. We introduce KORE-Map 1.0 (Korean medicine Omics Resource Extension Map), the first comprehensive resource of drug-derived transcriptome data for representative tonifying herbal medicines, effective in enhancing the immune system. This dataset aims to provide novel insights into the combinatorial mechanisms of these herbal medicines and to aid in the discovery of new therapeutic targets and indications for various incurable diseases.

Subject terms

High-throughput screening
RNA
https://doi.org/10.13039/501100003718 Korea Institute of Oriental Medicine (KIOM) KSN1722122 KSN1722122 KSN1722122 KSN1722122 KSN1722122 KSN1722122 KSN1722122 KSN1722122 Park Musun Park Sang-Min Lee Haeseung Kim Aeyung Yi Jin-Mu Cha Seongwon issue-copyright-statement© Springer Nature Limited 2024
==== Body
pmcBackground & Summary

Herbal medicine has long been used effectively for disease treatment in East Asia, notably in Korea and China1. These medicines exemplify combination therapy, employing multiple compounds to treat diseases2. Such therapy is effective in multifactorial diseases as it addresses multiple targets3, and offers the advantage of fewer side effects such as drug resistance4. Tonifying herbal medicine (THM), a form of combination therapy, not only targets diseases directly but also aims to activate immunity5. Thus, THMs have the capacity to: (1) treat multifactorial diseases, (2) minimize drug resistance, and (3) manage conditions challenging to treat with conventional drugs by bolstering the immune system. This study therefore proposes a new paradigm in disease treatment.

Despite its numerous benefits, research on THM faces limitations, particularly in accurately identifying the mechanism of action (MOA) due to the involvement of multiple compounds. One solution is drug-induced transcriptome analysis6,7, a method of measuring relative expression levels of mRNA after treating a cell line with a specific drug. This method reveals the MOA of drugs through complex pharmacodynamic processes as reflected in mRNA expression8, providing insights into the therapeutic effects of both single drugs and compound combinations. Furthermore, it facilitates the identification of transcriptomic signatures responsive to these treatments, confirming the therapeutic mechanism of each approach. Therefore, transcriptome analysis, by elucidating the MOA of combination therapies, is pivotal in revealing the therapeutic mechanisms of herbal medicine prescriptions beyond individual herbs.

To uncover the therapeutic MOA of herbal medicines using transcriptome analysis, it is important to consider a variety of factors. One such variable to consider is the choice of cell line9. As the dominantly expressed genes vary across cell lines, analysis across diverse cell lines essential for confirming the MOA of herbal medicines on multiple targets. Therefore, generating transcriptome data from various cell lines is an efficient strategy that not only validates the therapeutic effects of herbal medicines, but also elucidates their mechanisms of action in various organs. Another critical factor is drug concentration, which can significantly influence the MOA of drugs10. For example, Finasteride, a well-known treatment for benign prostatic hyperplasia, demonstrates different effects at varying doses—acting on prostate cancer at higher doses (5 mg/day), and delaying hair loss at lower doses (1 mg/day)11. This variability demonstrates the importance of dose determination in confirming the MOA of herbal medicines. Additionally, the choice of solvent plays a pivotal role7, as natural products comprise both hydrophilic and hydrophobic compounds. Therefore, the fractions of compounds extracted from water and ethanol solvents vary, indicating different processing mechanisms. To reveal the therapeutic mechanisms of herbal medicine, multiple variables must be considered. Transcriptome data, produced in consideration of these variables, are important in elucidating the MOA of the drug.

In this Data Descriptor, we introduce KORE-Map 1.0 (Korean Medicine Omics Resource Extension-Map), featuring THM-derived transcriptome data available on the NCBI gene expression omnibus (GEO) platform. The data were generated from four THMs commonly used in clinical practice, along with the 10 herbs constituting them. The transcript expression information was derived from both water and ethanol extracts of THMs and herbs, prepared at three different concentrations and applied to four representative human-derived cell lines (A549, HepG2, HT29, and SW1783). Utilizing the MGIEasy RNA directional library prep kit and the MGISEQ-2000 sequencing system, both widely used worldwide, not only facilitates easy data reuse but also ensures excellent compatibility with other transcriptome datasets. The THM-derived transcriptome data produced in this study could serve multiple purposes, such as aiding in the identification of therapeutic MOAs involving multiple compounds, which is crucial for understanding therapeutic mechanisms in multifactorial or comorbid conditions. In addition, the transcriptome data, spanning various cell lines and concentrations, hold potential for applications in drug repositioning, side effect detection, and more, by enabling the simultaneous evaluation of the effects of multiple compounds on multiple targets and organs.

Methods

Preparation of herbs

Dried medicinal plants, conforming to the Korean Pharmacopoeia standards, were provided by Kwangmyung-dang Medicinal Herbs Co., located in Ulsan, Republic of Korea. These samples underwent an organoleptic examination by Dr. Choi Goya, a herbal medicine organoleptic examination expert appointed by the Korea Food and Drug Administration. The identification to species level was accomplished through DNA barcode region sequencing. Voucher specimens have been stored at the Korean Herbarium of Standard Herbal Resources, within the Herbal Medicine Resources Research Center, at the Korea Institute of Oriental Medicine in Naju, Republic of Korea (Table 1). All herbs and extracts, which were sourced from the Oriental Medicine Resources Research Center (KIOM), are available online at https://oasis.kiom.re.kr/herblib.Table 1 Herbal Medicine Information and Yields.

Code Num.	Herbs	Extracts	
Medicinal name	Scientific Name	Family	Voucher number	D, Yield (%)	Voucher number	E, Yield (%)	Voucher number	
22-01D &E	Poria Sclerotium	Macrohyporia cocos (Schwein.) I.Johans. & Ryvarden	Polyporaceae	2-22-0221	2.5	3-22-0019	4.9	3-22-0020	
22-02D &E	Atractylodis Rhizoma Alba	Atractylodes macrocephala Koidz.	Compositae	2-22-0212	60.2	3-22-0017	19.7	3-22-0018	
22-03D &E	Glycyrrhizae Radix et Rhizoma	Glycyrrhiza glabra L.	Leguminosae	2-22-0211	28.8	3-22-0015	32.8	3-22-0016	
22-04D &E	Ginseng Radix	Panax ginseng C.A.Mey.	Araliaceae	2-22-0215	25.4	3-22-0013	20.9	3-22-0014	
22-05D &E	Angelica Root	Angelica acutiloba (Siebold & Zucc.) Kitag.	Umbelliferae	2-22-0217	46.2	3-22-0011	37.8	3-22-0012	
22-06D &E	Paeoniae Radix	Paeonia lactiflora Pall.	Paeoniaceae	2-22-0216	30.0	3-22-0009	28.7	3-22-0010	
22-07D &E	Cnidii Rhizoma	Ligusticum officinale (Makino) Kitag.	Umbelliferae	2-22-0218	42.2	3-22-0007	40.8	3-22-0008	
22-08D &E	Rehmanniae Radix Preparata	Rehmannia glutinosa (Gaertn.) DC.	Scrophulariaceae	2-22-0219	63.5	3-22-0004	67.1	3-22-0005	
22-09D &E	Cinnamomi Cortex	Cinnamomum cassia (L.) J.Presl	Lauraceae	2-22-0214	5.4	3-22-0003	11.5	3-22-0004	
22-10D &E	Astragali Radix	Astragalus mongholicus Bunge	Leguminosae	2-22-0213	19.6	3-22-0001	37.1	3-22-0002	
D: hot water extract; E: 70% ethanol extract.

Preparation method of hot water and 70% ethanol extracts of herbs and THMs

Hot water and 70% ethanol extracts of each plant were prepared and supplied by KOC Biotech Co., located in Daejeon, Republic of Korea. Initially, dried plants (1,000 g) were pulverized and extracted in 10 L of hot distilled water for 3 h using a reflux extraction system (MS-DM609; MTOPS, Seoul, Republic of Korea), or in 10 L of 70% ethanol for 1 h using an ultrasonication system (VCP-20, Lab companion, Dajeon, Republic of Korea) twice. The resulting extract solutions were filtered through a 5 µm cartridge filter, concentrated using a rotary evaporator (Ev-1020, SciLab, Seoul, Republic of Korea), and finally lyophilized in a freeze dryer (LP-20, Ilshin-Bio-Base, Dongducheon, Republic of Korea) to produce the final extracts. These extracts were then finely homogenized and packaged in glass bottles with desiccant silica gel. THMs were prepared by blending and homogenizing these extracts in accordance with the composition ratios and extract yields of the individual medicinal herbs, according to the Korean Pharmacopoeia (Table 2). For in vitro applications, extracts (100 mg) were vigorously vortexed for 30 min in 10 mL of phosphate-buffered saline (PBS; Thermo Fisher Scientific, Rockford, IL, USA) containing 2% DMSO. This mixture was then sterilized by filtration through a 0.22 µm membrane to obtain a stock solution (10 mg/mL), which was divided into small aliquots and stored at −80 °C until their use.Table 2 Tonifying Herbal Medicine Mixing Information.

THM	Herbs	Mass (g)	Proportion (%)	
Sagunja-tang	Poria Sclerotium	1.56	25	
Atractylodis Rhizoma Alba	1.56	25	
Glycyrrhizae Radix et Rhizoma	1.56	25	
Ginseng Radix	1.56	25	
Total mass	6.24	100	
Samul-tang	Angelica Root	1.33	25	
Paeoniae Radix	1.33	25	
Cnidii Rhizoma	1.33	25	
Rehmanniae Radix Preparata	1.33	25	
Total mass	5.32	100	
Palmul-tang	Glycyrrhizae Radix et Rhizoma	1.5	12.5	
Angelica Root	1.5	12.5	
Paeoniae Radix	1.5	12.5	
Atractylodis Rhizoma Alba	1.5	12.5	
Poria Sclerotium	1.5	12.5	
Rehmanniae Radix Preparata	1.5	12.5	
Ginseng Radix	1.5	12.5	
Cnidii Rhizoma	1.5	12.5	
Total mass	12	100	
Sibjeondeabo-tang	Ginseng Radix	1	10.53	
Atractylodis Rhizoma Alba	1	10.53	
Poria Sclerotium	1	10.53	
Glycyrrhizae Radix et Rhizoma	0.5	5.26	
Angelica Root	1	10.53	
Cnidii Rhizoma	1	10.53	
Paeoniae Radix	1	10.53	
Rehmanniae Radix Preparata	1	10.53	
Astragali Radix	1	10.53	
Cinnamomi Cortex	1	10.53	
Total mass	9.5	100.00	

Cell culture

All cell lines were purchased from the American Type Culture Collection (Manassas, VA, USA) and were cultured in a basal medium enriched with 10% heat-inactivated fetal bovine serum, 100 IU/mL penicillin, and 100 µg/mL streptomycin, all within a humidified incubator (Table 3). Cell confluence levels between 80–90% prompted the replacement of the growth medium every 3–4 days to maintain optimal growth conditions. To ensure the absence of mycoplasma contamination, the MycoAlert PLUS mycoplasma detection kit (Lonza, Rockland, ME, USA) was employed for regular testing.Table 3 Information of Cell Lines.

Cell Line	Information of cell line	
Growth Type	Description (Human)	Tissue	ATCC No.	Passage	Mycoplasma test	
A549	Adherent	Non-small cell lung carcinoma	Lung	CCL-185	85	negative	
HepG2	Adherent	Hepatocellular carcinoma cell line	Liver	HB-8065	81	negative	
HT29	Adherent	Colorectal adenocarcinoma	Colon	HTB-38	135	negative	
SW1783	Adherent	Astrocytoma	Brain	HTB-13	30	negative	
Cell Line	Culture condition	Cell density for test	
Medium for culture	Incubator	Maintanance	96 well for IC20	6 well for RNA samp.	
A549	RPMI + 10% FBS + antibiotics(P/S)	5% CO2, 37 °C	Trypsinize, 1:6, 3days	1 × 104cell/100uL	3 × 105 cells/3 mL	
HepG2	DMEM + 10% FBS + antibiotics(P/S)	5% CO2, 37 °C	Trypsinize, 1:4, 3days	2 × 104cell/100uL	6 × 105 cells/3 mL	
HT29	DMEM + 10% FBS + antibiotics(P/S)	5% CO2, 37 °C	Trypsinize, 1:4, 3days	2 × 104cell/100uL	6 × 105 cells/3 mL	
SW1783	Leibovitz’s L-15 + 10%FBS + antibiotics(P/S)	0% CO2, 37 °C	Trypsinize, 1:4, 4days	5 × 103cells/100uL	1.5 × 105cells/3 mL	

Drug treatment and total RNA preparation for RNA sequencing (RNA-seq) analysis

To determine the appropriate treatment drug concentrations, we performed cell cytotoxicity tests to investigate drug doses that maintained 80% cell viability (IC20s), which were then adopted as the maximal doses for RNA-seq data collection. For drugs whose IC20s could not be determined, the highest treatment concentrations were capped at 500 µg/mL for extracts, considering both their solubility and relevance for clinical application. To confirm the influence of concentration, cells were treated with three different concentrations using 1/5 serial dilutions, thereby exposing them to low, medium, and high doses. Positive control drugs such as wortmannin (Sigma, W1628), LY294002 (Sigma, L9908), and Thioridazine (Sigma, T9025) were incorporated into the assay for comparative analysis against the connectivity map (CMap) data. Cells treated with a 2% DMSO/PBS solution served as the vehicle control. One day before drug administration, cells were seeded into 6-well culture plates with 3 mL of growth medium. Following a 24 h treatment period, the cells were washed with ice-cold PBS, and total RNA was isolated using QIAzol RNA isolation reagents (Thermo Fisher Scientific) in accordance with the manufacturer’s instructions.

RNA-seq data generation and preprocessing

Total RNA (over 500 ng) from each sample was processed for the mRNA sequencing library using the MGIEasy RNA directional library prep kit (MGI Tech Co., Ltd., China), following the manufacturer’s instructions. The library concentration was quantified using the QuantiFluor® ssDNA System (Promega Corporation, WI, USA). The prepared DNA nanoball was sequenced on an MGISEQ system (MGI Tech Co., Ltd., China) employing 100 bp paired-end reads. The RNA-seq data quality was assessed using FastQC (v0.11.9). To remove common MGISEQ adapter sequences, TrimGalore (v0.6.6) was utilized. Trimmed reads were then mapped to the human reference genome assembly GRCh38 (hg38) using the STAR aligner (v2.7.3a) with default settings12. Gene transcript abundance, including expected read counts and transcripts per million, was quantified using RSEM (v1.3.3), with the gene annotation GRCh38.8413. The raw sequence data (FASTQ files) and the preprocessed expression values for each gene have been deposited in the GEO under accession numbers GSE244687, GSE244707, GSE244694, and GSE245912.

Differential gene expression analysis

Using the gene symbols of protein-coding genes, we utilized the collapseRows function from the WGCNA package (v.1.72-1)14, specifically designed to merge expression data for genes represented by multiple probes. This approach effectively reduces redundancy and potential noise, enhancing the clarity of subsequent analyses. Additionally, the filterByExpr function from the genefilter package (v.1.78.0)15 was utilized to exclude genes that failed to meet predetermined expression criteria across samples. This filtering ensured that only genes most likely to provide reliable and relevant signals were retained for analysis.

For evaluation of each set of treatment conditions—encompassing four cell lines, 14 herbs and THMs, two extraction methods, and three concentration levels— we conducted differential gene expression (DGE) analysis against the corresponding control samples. This analysis was performed using the Wald test statistic as implemented in the DESeq. 2 package (v.1.36.0)16. Differentially expressed genes (DEGs) were determined based on a fold-change threshold of 1.5 and an adjusted P-value of less than 0.05.

Clustering analysis

The fold-change values derived from the DGE analysis across all treatment conditions were clustered using the t-distributed stochastic neighbor embedding (t-SNE) algorithm. This machine learning technique, designed for dimensionality reduction, excels in visualizing high-dimensional datasets, making it a valuable tool for interpreting complex gene expression patterns. The analysis was conducted utilizing the Rtsne package (v.0.16), with the perplexity parameter set to 1017.

Comparisons with connectivity map transcriptome data

Connectivity Map data were obtained from the Clue.io platform(clue.io/data/CMap2020#LINCS2020). For our analysis, we selected level five gene expression signatures with high reproducibility, defined by moderated z-scores that met specific criteria (distil_cc_q75 > 0.5 and pct_self_rank_q25 > 0.05), to compare with our RNA-seq data. The R package CMapR (v1.8.0) was used to manipulate the level 5 GCTX file (level5_beta_trt_cp_ n720216 × 12328.gctx). Given the variance in gene expression profiling between our RNA-seq data and L1000 assays9 used in CMap, a direct comparison between gene expression values was difficult due to distinct distributions of expression values. To navigate this, we employed gene set enrichment analysis (GSEA) as an alternative method to explore the genome-wide perturbing effects of treatments such as wortmannin at the pathway level18. We utilized 2,229 gene sets from several databases—Hallmark, Biocarta, KEGG, REACTOME, PID, and Wikipathways—available through MSigDB (https://www.gsea-msigdb.org/gsea/msigdb). The analysis involved performing GSEA on all genes, ranked according to their Wald test statistics or level5 z-scores. To obtain the MSigDB gene sets and conduct GSEA, we utilized the R package MSigDBR (v7.5.1) and FGSEA (v3.18). From the GSEA results, we defined pathway activity score (PAC) as ‘sign (enrichment score) × -log10(p-value)’ value to quantify the significance level. PAC vectors of equal lengths (n = 2,229) were generated for both our dataset and the CMap dataset. Subsequently, we determined the Pearson correlation coefficient to assess the relationship between the PACs from our samples and those from CMap (Fig. 1).Fig. 1 Introduction to transcriptomic data production protocols and herbal drugs used. (A) Overview of standard operating procedures for producing standardized transcriptome data. (B) List of herbal drugs processed for transcriptome data production.

Data Records

All raw transcriptome data were uploaded to GEO in the FASTQ format using paired-end sequencing files. Each data file was presented in two fq.gz formats. Essential details such as the production method, adopted cell line, and dosage information were included in the metadata accompanying the GEO submission19–22 (Table 4). The dataset submitted to GEO comprised 1,092 RNA sequencing samples across 21 batches (Table 5, Supplementary Tables 1–6). Transcript samples were derived from four distinct cell lines: 270 (A549; accession number, GSE24468719), 270 (HepG2; accession number, GSE24468720), 273 (HT29; accession number, GSE24468721), and 279 (SW1783; accession number, GSE24468722). The difference in data volume between HT29 andSW1783 cell lines can be attributed to two factors:(1) A discrepancy in drugs used as positive controls, and (2) variations in the number of transcript production batches due to technical issues. Wortmannin, known for its anti-inflammatory properties, served as a universal positive control across all cell lines. Further, LY294002 produced in HT29 cells, and LY294002 and Thioridazine produced in SW1783 cells served as additional positive controls, contributing three and six samples, respectively. Consequently, six batches were specifically allocated for the SW1783 cell line, with an inclusion of three extra samples.Table 4 GEO accession number.

Cell lines	Accession number	URL	
A549	GSE24468719	https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE244687	
HepG2	GSE24470720	https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE244707	
HT29	GSE24469421	https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE244694	
SW1783	GSE24591222	https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE245912	

Table 5 Number of transcriptome data samples and batches per cell line.

Cell lines	#Samples	Batch	
A549	270	Batch A	60 samples (Vehicle, Wortmannin, 22-01D ~ 22-06D)	
Batch B	60 samples (Vehicle, 22-07D ~ 22-13D)	
Batch C	60 samples (Vehicle, 22-13D ~ 22-14D, 22-01E ~ 22-06E)	
Batch D	60 samples (Vehicle, 22-05E ~ 22-12E)	
Batch E	30 samples (Vehicle, 22-11E ~ 22-14E)	
HepG2	270	Batch A	60 samples (Vehicle, Wortmannin, 22-01D ~ 22-06D)	
Batch B	60 samples (Vehicle, 22-07D ~ 22-13D)	
Batch C	60 samples (Vehicle, 22-13D ~ 22-14D, 22-01E ~ 22-06E)	
Batch D	60 samples (Vehicle, 22-05E ~ 22-12E)	
Batch E	30 samples (Vehicle, 22-11E ~ 22-14E)	
HT29	273	Batch A	60 samples (Vehicle, Wortmannin, LY294002, 22-01D ~ 22-06D)	
Batch B	60 samples (Vehicle, 22-06D ~ 22-12D)	
Batch C	60 samples (Vehicle, 22-13D ~ 22-14D, 22-01E ~ 22-05E)	
Batch D	60 samples (Vehicle, 22-05E ~ 22-12E)	
Batch E	33 samples (Vehicle, 22-11E ~ 22-14E)	
SW1783	279	Batch A	60 samples (Vehicle, Wortmannin, LY294002, Thioridazine, 22-01D ~ 22-06D)	
Batch B	60 samples (Vehicle, 22-05D ~ 22-12D)	
Batch C	30 samples (Vehicle, 22-13D ~ 22-14D, 22-01E ~ 22-02E)	
Batch D	9 samples (Vehicle, 22-01E ~ 22-02E)	
Batch E	60 samples (Vehicle, 22-01E ~ 22-08E)	
Batch F	60 samples (Vehicle, 22-08E ~ 22-14E)	

Technical Validation

RNA quality and integrity

To ensure the suitability of samples for downstream sequencing, RNA quality and integrity were thoroughly evaluated. The optical density at 260 and 280 nm was measured using the Trinean Dropsense™96 micro-volume reader. The A260/A280 ratio serves as an estimate of RNA purity, with values exceeding 1.8 indicating relatively high purity. Our analysis revealed that the RNA samples typically exhibited a ratio close to 1.8, signifying a substantial proportion of pure RNA (Fig. 2a). Furthermore, the 28S/18S rRNA ratio and the RNA integrity number (RIN) were measured using an agilent bioanalyzer DNA chip to assess the extent of RNA degradation. All RNA samples demonstrated a 28S/18S ratio approximately equal to 2 and an RIN value of 7 or above, reflecting high RNA quality and integrity (Fig. 2b). These results suggest that the RNA is of suitable quality for downstream RNA sequencing23.Fig. 2 Quality assessment of RNA samples. (A) The A260/A280 ratio for individual samples grouped by four cell lines. The minimum of values widely recognized as indicative of high purity RNA is represented by the dotted line. (B) The 28 s/18 s rRNA ratio (left) and RNA integrity number (right) for individual samples grouped by four cell lines. Each dotted line represents a minimum value widely known to reflect high RNA quality and integrity.

Quality of RNA-seq data

The quality of the raw RNA-seq data was assessed using FastQC (v0.11.9), a software that generates a detailed report, including metrics such as per-base quality scores and GC content distribution. A representative FastQC report indicated that the overall read quality was high (Fig. 3b). Similar quality metrics were observed in all other FastQC reports, qualifying them for further analysis. To obtain clean data, adapter sequences and low-quality bases (Phred score below 20) were removed using TrimGalore (v0.6.6). As a result, a high percentage of reads, with a median of 96.67%, were successfully and uniquely mapped to the human reference genome GRCh38 (hg38) (Fig. 3b)24.Fig. 3 Quality evaluation of RNA sequencing data. (A) Representative FastQC report showing per sequence quality scores (left) and GC content (right) for A549 cell line treated with dimethyl sulfoxide (DMSO). (B) Summary of unmapped, multiple-mapped, and uniquely mapped reads against the reference genome for each cell line.

Biological and technical reproducibility

To ensure the reproducibility of our RNA-seq data, we quantified biological and technical batch effects by analyzing expression levels (TPM values for 19,826 protein-coding genes). Initially, the biological reproducibility was assessed through the analysis of three independent biological replicates for each treatment condition; cell line, treatment, and dose. Each replicate involved separate RNA extraction, RNA-seq library preparation, and sequencing processes. We calculated the pairwise Pearson’s correlation coefficient among replicates to quantify their similarities. This revealed a high degree of biological reproducibility, with an average correlation coefficient of 0.994 across all conditions. Furthermore, 97.8% of the conditions exhibited an average expression level correlation exceeding 0.95 across the three replicates (Fig. 4a).Fig. 4 Replicability of RNA-seq profiles. (A) Distribution of Pearson’s correlation coefficients for replicates (yellow) versus different samples (gray). (B) Heatmap of Pearson’s correlation coefficients among replicate samples across various sequencing batches.

Technical reproducibility was subsequently evaluated to address potential batch effects arising from sequencing. Since a single sequencing lane can accommodate up to 60 samples, we distributed samples from the same cell line across six different sequencing batches (A to F). Control samples treated with the vehicle were included and sequenced in all six batches, to assess batch effects. The correlation coefficients between control samples from different batches were calculated, indicating minimal batch effects. Notably, all control samples exhibited high correlation coefficients (>0.99) with samples sequenced in different batches (Fig. 4b).

Comparisons with external drug-induced transcriptome data

To verify the reliability and reproducibility of our RNA-seq data, we compared our drug-induced transcriptome profiles to those generated by the CMap 9, a comprehensive database featuring gene expression profiles of human cell lines treated with various bioactive compounds. We chose wortmannin, an established positive control that is also included in the CMap dataset, as a benchmark for our analyses across three cell lines: A549, HEPG2, and HT29.

To facilitate direct comparison between transcriptome data generated from different platforms, we converted gene-level expression values to pathway-level scores. This approach aggregates the expression changes across genes within 2,229 well-defined biological pathways, providing a more robust and interpretable measure of pathway activation or inhibition.

We then compared the pathway-level scores resulting from wortmannin treatment in our study with those generated by CMap. The pathway-level scores from our wortmannin treatment analysis showed a high correlation with those obtained from CMap (Fig. 5). This notable concordance serves as strong evidence of the reliability and reproducibility of our RNA-seq data, affirming its ability to capture drug-induced changes in cellular pathways.Fig. 5 Transcriptome data comparisons with CMap. Distribution of Pearson’s correlation coefficients for samples under the same treatment condition (yellow) versus different conditions (gray).

Supplementary information

Supplementary Table 1

Supplementary Table 2

Supplementary Table 3

Supplementary Table 4

Supplementary Table 5

Supplementary Table 6

Supplementary information

The online version contains supplementary material available at 10.1038/s41597-024-03734-x.

Acknowledgements

Total RNA isolation and RNA sequencing were conducted by LAS Co. Ltd., South Korea. This study was supported by Grant number KSN1722122 from the Korea Institute of Oriental Medicine.

Author contributions

All authors contributed to the study. M.P. and J.Y. were responsible for data generation. M.P., S.P., and H.L. analyzed the data. M.P. uploaded and curated data. S.P. and H.L. visualized the data and performed formal analysis. M.P., S.P., H.L., J.Y. and S.C. drafted the manuscript. All authors reviewed and edited the manuscript.

Code availability

All software used to analyze the RNA-seq data, along with their parameters are clearly described in the Methods section. Unless specified otherwise, default parameter settings recommended by the developers were used. The curation and validation of the dataset were conducted using custom R code, as detailed in the Materials and Methods section. The source code specifically employed for preprocessing the RNA-seq data is publicly available through the GitHub repository (https://github.com/LeeLab-Sysbio/KOREMAP.v1). Researchers are encouraged to cite this paper when utilizing the RNA-seq data uploaded in GEO.

Competing interests

The authors declare no competing interests.

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

These authors contributed equally: Musun Park, Sang-Min Park, Haeseung Lee.
==== Refs
References

1. Tyler VE Herbal medicine: from the past to the future Public Health Nutr 2000 3 447 452 10.1017/s1368980000000525 11276292
Tyler, V. E. Herbal medicine: from the past to the future. Public Health Nutr 3, 447–452, 10.1017/s1368980000000525 (2000).11276292 10.1017/s1368980000000525
2. Gouws C Steyn D Du Plessis L Steenekamp J Hamman JH Combination therapy of Western drugs and herbal medicines: recent advances in understanding interactions involving metabolism and efflux Expert Opin Drug Metab Toxicol 2012 8 973 984 10.1517/17425255.2012.691966 22612723
Gouws, C., Steyn, D., Du Plessis, L., Steenekamp, J. & Hamman, J. H. Combination therapy of Western drugs and herbal medicines: recent advances in understanding interactions involving metabolism and efflux. Expert Opin Drug Metab Toxicol 8, 973–984, 10.1517/17425255.2012.691966 (2012).22612723 10.1517/17425255.2012.691966
3. Rosini M Polypharmacology: the rise of multitarget drugs over combination therapies Future Med Chem 2014 6 485 487 10.4155/fmc.14.25 24649950
Rosini, M. Polypharmacology: the rise of multitarget drugs over combination therapies. Future Med Chem 6, 485–487, 10.4155/fmc.14.25 (2014).24649950 10.4155/fmc.14.25
4. Basanta D Gatenby RA Anderson AR Exploiting evolution to treat drug resistance: combination therapy and the double bind Mol Pharm 2012 9 914 921 10.1021/mp200458e 22369188
Basanta, D., Gatenby, R. A. & Anderson, A. R. Exploiting evolution to treat drug resistance: combination therapy and the double bind. Mol Pharm 9, 914–921, 10.1021/mp200458e (2012).22369188 10.1021/mp200458e
5. Leung HY Ko KM Differential Effects of Yin-and Yang-Chinese Tonifying Herbs on Innate and Adaptive Immunity Chinese Medicine 2023 14 68 78 10.4236/cm.2023.142004
Leung, H. Y. & Ko, K. M. Differential Effects of Yin-and Yang-Chinese Tonifying Herbs on Innate and Adaptive Immunity. Chinese Medicine 14, 68–78 (2023).10.4236/cm.2023.142004
6. Baek SJ Identification of a novel anticancer mechanism of Paeoniae Radix extracts based on systematic transcriptome analysis Biomed Pharmacother 2022 148 112748 10.1016/j.biopha.2022.112748 35219117
Baek, S. J. et al. Identification of a novel anticancer mechanism of Paeoniae Radix extracts based on systematic transcriptome analysis. Biomed Pharmacother 148, 112748, 10.1016/j.biopha.2022.112748 (2022).35219117 10.1016/j.biopha.2022.112748
7. Park SM Systematic transcriptome analysis reveals molecular mechanisms and indications of bupleuri radix Front Pharmacol 2022 13 1010520 10.3389/fphar.2022.1010520 36304143
Park, S. M. et al. Systematic transcriptome analysis reveals molecular mechanisms and indications of bupleuri radix. Front Pharmacol 13, 1010520, 10.3389/fphar.2022.1010520 (2022).36304143 10.3389/fphar.2022.1010520
8. Lee H Kang S Kim W Drug Repositioning for Cancer Therapy Based on Large-Scale Drug-Induced Transcriptional Signatures PLoS One 2016 11 e0150460 10.1371/journal.pone.0150460 26954019
Lee, H., Kang, S. & Kim, W. Drug Repositioning for Cancer Therapy Based on Large-Scale Drug-Induced Transcriptional Signatures. PLoS One 11, e0150460, 10.1371/journal.pone.0150460 (2016).26954019 10.1371/journal.pone.0150460
9. Subramanian A A Next Generation Connectivity Map: L1000 Platform and the First 1,000,000 Profiles Cell 2017 171 1437 1452 e1417 10.1016/j.cell.2017.10.049 29195078
Subramanian, A. et al. A Next Generation Connectivity Map: L1000 Platform and the First 1,000,000 Profiles. Cell 171, 1437–1452 e1417, 10.1016/j.cell.2017.10.049 (2017).29195078 10.1016/j.cell.2017.10.049
10. Waldmann T Design principles of concentration-dependent transcriptome deviations in drug-exposed differentiating stem cells Chem Res Toxicol 2014 27 408 420 10.1021/tx400402j 24383497
Waldmann, T. et al. Design principles of concentration-dependent transcriptome deviations in drug-exposed differentiating stem cells. Chem Res Toxicol 27, 408–420, 10.1021/tx400402j (2014).24383497 10.1021/tx400402j
11. Libecco JF Bergfeldt WF Finasteride in the treatment of alopecia Expert Opin Pharmaco 2004 5 933 940 10.1517/14656566.5.4.933
Libecco, J. F. & Bergfeldt, W. F. Finasteride in the treatment of alopecia. Expert Opin Pharmaco 5, 933–940, 10.1517/14656566.5.4.933 (2004).10.1517/14656566.5.4.933
12. Dobin A STAR: ultrafast universal RNA-seq aligner Bioinformatics 2013 29 15 21 10.1093/bioinformatics/bts635 23104886
Dobin, A. et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21, 10.1093/bioinformatics/bts635 (2013).23104886 10.1093/bioinformatics/bts635
13. Li, B. & Dewey, C. N. RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. Bmc Bioinformatics 1210.1186/1471-2105-12-323 (2011).
14. Langfelder, P. & Horvath, S. WGCNA: an R package for weighted correlation network analysis. Bmc Bioinformatics 910.1186/1471-2105-9-559 (2008).
15. Gentleman, R., Carey, V., Huber, W. & Hahne, F. Genefilter: methods for filtering genes from high-throughput experiments. R package version 1 (2015).
16. Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq. 2. Genome Biol 1510.1186/s13059-014-0550-8 (2014).
17. van der Maaten L Hinton G Visualizing Data using t-SNE J Mach Learn Res 2008 9 2579 2605
van der Maaten, L. & Hinton, G. Visualizing Data using t-SNE. J Mach Learn Res 9, 2579–2605 (2008).
18. Subramanian A Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles Proc Natl Acad Sci USA 2005 102 15545 15550 10.1073/pnas.0506580102 16199517
Subramanian, A. et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci USA 102, 15545–15550, 10.1073/pnas.0506580102 (2005).16199517 10.1073/pnas.0506580102
19. Park M Cha S Yi J-M 2023 Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_A549 GSE244687
Park, M., Cha, S. & Yi, J.-M. Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_A549, https://identifiers.org/GEO:GSE244687 (2023).
20. Park M Cha S Yi J-M 2023 Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_HepG2 GSE244707
Park, M., Cha, S. & Yi, J.-M. Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_HepG2, https://identifiers.org/GEO:GSE244707 (2023).
21. Park M Cha S Yi J-M 2023 Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_HT29 GSE244694
Park, M., Cha, S. & Yi, J.-M. Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_HT29, https://identifiers.org/GEO:GSE244694 (2023).
22. Park M Cha S Yi J-M 2023 Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_SW1783 GSE245912
Park, M., Cha, S. & Yi, J.-M. Digital transformation of herbal medicine: Conversion to biological entity data using tonifying herbal medicine-induced transcriptome sequencing_Tonifying_SW1783, https://identifiers.org/GEO:GSE245912 (2023).
23. Park M 2024 KORE-map_Sample quality 10.6084/m9.figshare.26019445.v1
Park, M. et al. KORE-map_Sample quality, 10.6084/m9.figshare.26019445.v1 (2024).10.6084/m9.figshare.26019445.v1
24. Park M 2024 KORE-map_Sequencing quality 10.6084/m9.figshare.26020525.v1
Park, M. et al. KORE-map_Sequencing quality, 10.6084/m9.figshare.26020525.v1 (2024).10.6084/m9.figshare.26020525.v1
