ATAC-seq with unique molecular identifiers improves quantification and footprinting

Time-course transcriptome and chromatin accessibility analyses reveal the dynamic transcriptional regulation shaping spikelet hull size

The grains of rice (Oryza sativa) are enclosed by a spikelet hull comprising the lemma and palea. Development of the spikelet hull determines the storage capacity of the grain, thus affecting grain yield and quality. Although multiple signaling pathways controlling grain size have been identified, the transcriptional regulatory mechanisms underlying grain development remain limited. Here, we used RNA-seq and ATAC-seq to characterize the transcription and chromatin accessibility dynamics during the development of spikelet hulls. A time-course analysis showed that more than half of the genes were sequentially expressed during hull development and that the accessibility of most open chromatin regions (OCRs) changed moderately, although some regions positively or negatively affected the expression of their closest genes. We revealed a crucial role of GROWTH-REGULATING FACTORs in shaping grain size by influencing multiple metabolic and signaling pathways, and a coordinated transcriptional regulation in response to auxin and cytokinin signaling. We also demonstrated the function of SCL6-IIb, a member of the GRAS family transcription factors, in regulating grain size, with SCL6-IIb expression being activated by SQUAMOSA PROMOTER BINDING PROTEIN-LIKE 18 (OsSPL18). When we edited the DNA sequences within OCRs upstream of the start codon of BRASSINAZOLE-RESISTANT 1 (BZR1) and SCL6-IIb, we generated multiple mutant lines with longer grains. These findings offer a comprehensive overview of the cis-regulatory landscape involved in forming grain capacity and a valuable resource for exploring the regulatory network behind grain development.

Comprehensive mapping and modelling of the rice regulome landscape unveils the regulatory architecture underlying complex traits

Unraveling the regulatory mechanisms that govern complex traits is pivotal for advancing crop improvement. Here we present a comprehensive regulome atlas for rice (Oryza sativa), charting the chromatin accessibility across 23 distinct tissues from three representative varieties. Our study uncovers 117,176 unique open chromatin regions (OCRs), accounting for ~15% of the rice genome, a notably higher proportion compared to previous reports in plants. Integrating RNA-seq data from matched tissues, we confidently predict 59,075 OCR-to-gene links, with enhancers constituting 69.54% of these associations, including many known enhancer-to-gene links. Leveraging this resource, we re-evaluate genome-wide association study results and discover a previously unknown function of OsbZIP06 in seed germination, which we subsequently confirm through experimental validation. We optimize deep learning models to decode regulatory grammar, achieving robust modeling of tissue-specific chromatin accessibility. This approach allows to predict cross-variety regulatory dynamics from genomic sequences, shedding light on the genetic underpinnings of cis-regulatory divergence and morphological disparities between varieties. Overall, our study establishes a foundational resource for rice functional genomics and precision molecular breeding, providing valuable insights into regulatory mechanisms governing complex traits.

Characterization of the DNA accessibility of chloroplast genomes in grasses

Although the chloroplast genome (cpDNA) of higher plants is known to exist as a large protein-DNA complex called 'plastid nucleoid', researches on its DNA state and regulatory elements are limited. In this study, we performed the assay for transposase-accessible chromatin sequencing (ATAC-seq) on five common tissues across five grasses, and found that the accessibility of different regions in cpDNA varied widely, with the transcribed regions being highly accessible and accessibility patterns around gene start and end sites varying depending on the level of gene expression. Further analysis identified a total of 3970 putative protein binding footprints on cpDNAs of five grasses. These footprints were enriched in intergenic regions and co-localized with known functional elements. Footprints and their flanking accessibility varied dynamically among tissues. Cross-species analysis showed that footprints in coding regions tended to overlap non-degenerate sites and contain a high proportion of highly conserved sites, indicating that they are subject to evolutionary constraints. Taken together, our results suggest that the accessibility of cpDNA has biological implications and provide new insights into the transcriptional regulation of chloroplasts.

Principles of digital sequencing using unique molecular identifiers

Massively parallel sequencing technologies have long been used in both basic research and clinical routine. The recent introduction of digital sequencing has made previously challenging applications possible by significantly improving sensitivity and specificity to now allow detection of rare sequence variants, even at single molecule level. Digital sequencing utilizes unique molecular identifiers (UMIs) to minimize sequencing-induced errors and quantification biases. Here, we discuss the principles of UMIs and how they are used in digital sequencing. We outline the properties of different UMI types and the consequences of various UMI approaches in relation to experimental protocols and bioinformatics. Finally, we describe how digital sequencing can be applied in specific research fields, focusing on cancer management where it can be used in screening of asymptomatic individuals, diagnosis, treatment prediction, prognostication, monitoring treatment efficacy and early detection of treatment resistance as well as relapse.

Single-cell ATAC sequencing identifies sleepy macrophages during reciprocity of cytokines in L. major infection

The hallmark characteristic of macrophages lies in their inherent plasticity, allowing them to adapt to dynamic microenvironments. Leishmania strategically modulates the phenotypic plasticity of macrophages, creating a favorable environment for intracellular survival and persistent infection through regulatory cytokine such as interleukin (IL)-10. Nevertheless, these effector cells can counteract infection by modulating crucial cytokines like IL-12 and key components involved in its production. Using sophisticated tool of single-cell assay for transposase accessible chromatin (ATAC) sequencing, we systematically examined the regulatory axis of IL-10 and IL-12 in a time-dependent manner during Leishmania major infection in macrophages Our analysis revealed the cellular heterogeneity post-infection with the regulators of IL-10 and IL-12, unveiling a reciprocal relationship between these cytokines. Notably, our significant findings highlighted the presence of sleepy macrophages and their pivotal role in mediating reciprocity between IL-10 and IL-12. To summarize, the roles of cytokine expression, transcription factors, cell cycle, and epigenetics of host cell machinery were vital in identification of sleepy macrophages, which is a transient state where transcription factors controlled the epigenetic remodeling and expression of genes involved in pro-inflammatory cytokine expression and recruitment of immune cells.IMPORTANCELeishmaniasis is an endemic affecting 99 countries and territories globally, as outlined in the 2022 World Health Organization report. The disease's severity is compounded by compromised host immune systems, emphasizing the pivotal role of the interplay between parasite and host immune factors in disease regulation. In instances of cutaneous leishmaniasis induced by L. major, macrophages function as sentinel cells. Our findings indicate that the plasticity and phenotype of macrophages can be modulated to express a cytokine profile involving IL-10 and IL-12, mediated by the regulation of transcription factors and their target genes post-L. major infection in macrophages. Employing sophisticated methodologies such as single-cell ATAC sequencing and computational genomics, we have identified a distinctive subset of macrophages termed "sleepy macrophages." These macrophages exhibit downregulated housekeeping genes while expressing a unique set of variable features. This data set constitutes a valuable resource for comprehending the intricate host-parasite interplay during L. major infection.

A high-performance computational workflow to accelerate GATK SNP detection across a 25-genome dataset

BACKGROUND

Single-nucleotide polymorphisms (SNPs) are the most widely used form of molecular genetic variation studies. As reference genomes and resequencing data sets expand exponentially, tools must be in place to call SNPs at a similar pace. The genome analysis toolkit (GATK) is one of the most widely used SNP calling software tools publicly available, but unfortunately, high-performance computing versions of this tool have yet to become widely available and affordable.

RESULTS

Here we report an open-source high-performance computing genome variant calling workflow (HPC-GVCW) for GATK that can run on multiple computing platforms from supercomputers to desktop machines. We benchmarked HPC-GVCW on multiple crop species for performance and accuracy with comparable results with previously published reports (using GATK alone). Finally, we used HPC-GVCW in production mode to call SNPs on a "subpopulation aware" 16-genome rice reference panel with ~ 3000 resequenced rice accessions. The entire process took ~ 16 weeks and resulted in the identification of an average of 27.3 M SNPs/genome and the discovery of ~ 2.3 million novel SNPs that were not present in the flagship reference genome for rice (i.e., IRGSP RefSeq).

CONCLUSIONS

This study developed an open-source pipeline (HPC-GVCW) to run GATK on HPC platforms, which significantly improved the speed at which SNPs can be called. The workflow is widely applicable as demonstrated successfully for four major crop species with genomes ranging in size from 400 Mb to 2.4 Gb. Using HPC-GVCW in production mode to call SNPs on a 25 multi-crop-reference genome data set produced over 1.1 billion SNPs that were publicly released for functional and breeding studies. For rice, many novel SNPs were identified and were found to reside within genes and open chromatin regions that are predicted to have functional consequences. Combined, our results demonstrate the usefulness of combining a high-performance SNP calling architecture solution with a subpopulation-aware reference genome panel for rapid SNP discovery and public deployment.

Capture of regulatory factors via CRISPR-dCas9 for mechanistic analysis of fine-tuned SERRATE expression in Arabidopsis

SERRATE (SE) plays an important role in many biological processes and under biotic stress resistance. However, little about the control of SE has been clarified. Here we present a method named native chromatin-associated proteome affinity by CRISPR-dCas9 (CASPA-dCas9) to holistically capture native regulators of the SE locus. Several key regulatory factors including PHYTOCHROME RAPIDLY REGULATED 2 (PAR2), WRKY DNA-binding protein 19 (WRKY19) and the MYB-family protein MYB27 of SE are identified. MYB27 recruits the long non-coding RNA-PRC2 (SEAIR-PRC2) complex for H3K27me3 deposition on exon 1 of SE and subsequently represses SE expression, while PAR2-MYB27 interaction inhibits both the binding of MYB27 on the SE promoter and the recruitment of SEAIR-PRC2 by MYB27. The interaction between PAR2 and MYB27 fine-tunes the SE expression level at different developmental stages. In addition, PAR2 and WRKY19 synergistically promote SE expression for pathogen resistance. Collectively, our results demonstrate an efficient method to capture key regulators of target genes and uncover the precise regulatory mechanism for SE.

cisDynet: An integrated platform for modeling gene-regulatory dynamics and networks

Chromatin accessibility sequencing has been widely used for uncovering genetic regulatory mechanisms and inferring gene regulatory networks. However, effectively integrating large-scale chromatin accessibility datasets has posed a significant challenge. This is due to the lack of a comprehensive end-to-end solution, as many existing tools primarily emphasize data preprocessing and overlook downstream analyses. To bridge this gap, we have introduced cisDynet, a holistic solution that combines streamlined data preprocessing using Snakemake and R functions with advanced downstream analysis capabilities. cisDynet excels in conventional data analyses, encompassing peak statistics, peak annotation, differential analysis, motif enrichment analysis, and more. Additionally, it allows to perform sophisticated data exploration, such as tissue-specific peak identification, time course data modeling, integration of RNA-seq data to establish peak-to-gene associations, constructing regulatory networks, and conducting enrichment analysis of genome-wide association study (GWAS) variants. As a proof of concept, we applied cisDynet to reanalyze comprehensive ATAC-seq datasets across various tissues from the Encyclopedia of DNA Elements (ENCODE) project. The analysis successfully delineated tissue-specific open chromatin regions (OCRs), established connections between OCRs and target genes, and effectively linked these discoveries with 1861 GWAS variants. Furthermore, cisDynet was instrumental in dissecting the time course open chromatin data of mouse embryonic development, revealing the dynamic behavior of OCRs over developmental stages and identifying key transcription factors governing differentiation trajectories. In summary, cisDynet offers researchers a user-friendly solution that minimizes the need for extensive coding, ensures the reproducibility of results, and greatly simplifies the exploration of epigenomic data.

Integrated 3D genome, epigenome and transcriptome analyses reveal transcriptional coordination of circadian rhythm in rice

Photoperiods integrate with the circadian clock to coordinate gene expression rhythms and thus ensure plant fitness to the environment. Genome-wide characterization and comparison of rhythmic genes under different light conditions revealed delayed phase under constant darkness (DD) and reduced amplitude under constant light (LL) in rice. Interestingly, ChIP-seq and RNA-seq profiling of rhythmic genes exhibit synchronous circadian oscillation in H3K9ac modifications at their loci and long non-coding RNAs (lncRNAs) expression at proximal loci. To investigate how gene expression rhythm is regulated in rice, we profiled the open chromatin regions and transcription factor (TF) footprints by time-series ATAC-seq. Although open chromatin regions did not show circadian change, a significant number of TFs were identified to rhythmically associate with chromatin and drive gene expression in a time-dependent manner. Further transcriptional regulatory networks mapping uncovered significant correlation between core clock genes and transcription factors involved in light/temperature signaling. In situ Hi-C of ZT8-specific expressed genes displayed highly connected chromatin association at the same time, whereas this ZT8 chromatin connection network dissociates at ZT20, suggesting the circadian control of gene expression by dynamic spatial chromatin conformation. These findings together implicate the existence of a synchronization mechanism between circadian H3K9ac modifications, chromatin association of TF and gene expression, and provides insights into circadian dynamics of spatial chromatin conformation that associate with gene expression rhythms.

Protocol for Bulk-ATAC sequencing in head and neck squamous cell carcinoma

The transposase-accessible chromatin using sequencing (ATAC-seq) offers a simplified approach to detect chromatin changes in cancer cells after genetic intervention and drug treatment. Here, we present an optimized ATAC-seq protocol to elucidate chromatin accessibility changes at the epigenetic level in head and neck squamous cell carcinoma cells. We describe steps for cell lysate preparation, transposition, and tagmentation, followed by library amplification and purification. We then detail next-generation sequencing and data analysis. For complete details on the use and execution of this protocol, please refer to Buenrostro et al.,¹ Chen et al.,².

Systematic annotation of conservation states provides insights into regulatory regions in rice

Plant genomes contain a large fraction of noncoding sequences. The discovery and annotation of conserved noncoding sequences (CNSs) in plants is an ongoing challenge. Here we report the application of comparative genomics to systematically identify CNSs in 50 well-annotated Gramineae genomes using rice (Oryza sativa) as the reference. We conduct multiple-way whole-genome alignments to the rice genome. The rice genome is annotated as 20 conservation states (CSs) at single-nucleotide resolution using a multivariate hidden Markov model (ConsHMM) based on the multiple-genome alignments. Different states show distinct enrichments for various genomic features, and the conservation scores of CSs are highly correlated with the level of associated chromatin accessibility. We find that at least 33.5% of the rice genome is highly under selection, with more than 70% of the sequence lying outside of coding regions. A catalog of 855,366 regulatory CNSs is generated, and they significantly overlapped with putative active regulatory elements such as promoters, enhancers, and transcription factor binding sites. Collectively, our study provides a resource for elucidating functional noncoding regions of the rice genome and an evolutionary aspect of regulatory sequences in higher plants.

Genomic innovation and regulatory rewiring during evolution of the cotton genus Gossypium

Phenotypic diversity and evolutionary innovation ultimately trace to variation in genomic sequence and rewiring of regulatory networks. Here, we constructed a pan-genome of the Gossypium genus using ten representative diploid genomes. We document the genomic evolutionary history and the impact of lineage-specific transposon amplification on differential genome composition. The pan-3D genome reveals evolutionary connections between transposon-driven genome size variation and both higher-order chromatin structure reorganization and the rewiring of chromatin interactome. We linked changes in chromatin structures to phenotypic differences in cotton fiber and identified regulatory variations that decode the genetic basis of fiber length, the latter enabled by sequencing 1,005 transcriptomes during fiber development. We showcase how pan-genomic, pan-3D genomic and genetic regulatory data serve as a resource for delineating the evolutionary basis of spinnable cotton fiber. Our work provides insights into the evolution of genome organization and regulation and will inform cotton improvement by enabling regulome-based approaches.

Comprehensive transcriptional variability analysis reveals gene networks regulating seed oil content of Brassica napus

BACKGROUND

Regulation of gene expression plays an essential role in controlling the phenotypes of plants. Brassica napus (B. napus) is an important source for the vegetable oil in the world, and the seed oil content is an important trait of B. napus.

RESULTS

We perform a comprehensive analysis of the transcriptional variability in the seeds of B. napus at two developmental stages, 20 and 40 days after flowering (DAF). We detect 53,759 and 53,550 independent expression quantitative trait loci (eQTLs) for 79,605 and 76,713 expressed genes at 20 and 40 DAF, respectively. Among them, the local eQTLs are mapped to the adjacent genes more frequently. The adjacent gene pairs are regulated by local eQTLs with the same open chromatin state and show a stronger mode of expression piggybacking. Inter-subgenomic analysis indicates that there is a feedback regulation for the homoeologous gene pairs to maintain partial expression dosage. We also identify 141 eQTL hotspots and find that hotspot87-88 co-localizes with a QTL for the seed oil content. To further resolve the regulatory network of this eQTL hotspot, we construct the XGBoost model using 856 RNA-seq datasets and the Basenji model using 59 ATAC-seq datasets. Using these two models, we predict the mechanisms affecting the seed oil content regulated by hotspot87-88 and experimentally validate that the transcription factors, NAC13 and SCL31, positively regulate the seed oil content.

CONCLUSIONS

We comprehensively characterize the gene regulatory features in the seeds of B. napus and reveal the gene networks regulating the seed oil content of B. napus.

Mapping open chromatin by ATAC-seq in bread wheat

Gene transcription is largely regulated by cis-regulatory elements. Assay for Transposase-Accessible Chromatin using sequencing (ATAC-seq) is an emerging technology that can accurately map cis-regulatory elements in animals and plants. However, the presence of cell walls and chloroplasts in plants hinders the extraction of high-quality nuclei, thereby affects the quality of ATAC-seq data. Meanwhile, it is tricky to perform ATAC-seq with different tissue types, especially for those with limited size and amount. Moreover, with rapid growth of ATAC-seq datasets from plants, powerful and easy-to-use data analysis pipelines for ATAC-seq, especially for wheat is lacking. Here, we provided an all-in-one solution for mapping open chromatin in wheat including both experimental and data analysis procedure. We efficiently obtained nuclei with less cell debris from various wheat tissues. High-quality ATAC-seq data from young spike and ovary, which are hard to harvest were generated. We determined that the saturation sequencing depth of wheat ATAC-seq is about 16 Gb. Particularly, we developed a powerful and easy-to-use online pipeline to analyze the wheat ATAC-seq data and this pipeline can be easily extended to other plant species. The method developed here will facilitate plant regulatory genome study not only for wheat but also for other plant species.

An inferred functional impact map of genetic variants in rice

Interpreting the functional impacts of genetic variants (GVs) is an important challenge for functional genomic studies in crops and next-generation breeding. Previous studies in rice (Oryza sativa) have focused mainly on the identification of GVs, whereas systematic functional annotation of GVs has not yet been performed. Here, we present a functional impact map of GVs in rice. We curated haplotype information for 17 397 026 GVs from sequencing data of 4726 rice accessions. We quantitatively evaluated the effects of missense mutations in coding regions in each haplotype based on the conservation of amino acid residues and obtained the effects of 918 848 non-redundant missense GVs. Furthermore, we generated high-quality chromatin accessibility (CA) data from six representative rice tissues and used these data to train deep convolutional neural network models to predict the impacts of 5 067 405 GVs for CA in regulatory regions. We characterized the functional properties and tissue specificity of the GV effects and found that large-effect GVs in coding and regulatory regions may be subject to selection in different directions. Finally, we demonstrated how the functional impact map could be used to prioritize causal variants in mapping populations. This impact map will be a useful resource for accelerating gene cloning and functional studies in rice, and can be freely queried in RiceVarMap V2.0 (http://ricevarmap.ncpgr.cn).

PlantDeepSEA, a deep learning-based web service to predict the regulatory effects of genomic variants in plants

Characterizing regulatory effects of genomic variants in plants remains a challenge. Although several tools based on deep-learning models and large-scale chromatin-profiling data have been available to predict regulatory elements and variant effects, no dedicated tools or web services have been reported in plants. Here, we present PlantDeepSEA as a deep learning-based web service to predict regulatory effects of genomic variants in multiple tissues of six plant species (including four crops). PlantDeepSEA provides two main functions. One is called Variant Effector, which aims to predict the effects of sequence variants on chromatin accessibility. Another is Sequence Profiler, a utility that performs 'in silico saturated mutagenesis' analysis to discover high-impact sites (e.g., cis-regulatory elements) within a sequence. When validated on independent test sets, the area under receiver operating characteristic curve of deep learning models in PlantDeepSEA ranges from 0.93 to 0.99. We demonstrate the usability of the web service with two examples. PlantDeepSEA could help to prioritize regulatory causal variants and might improve our understanding of their mechanisms of action in different tissues in plants. PlantDeepSEA is available at http://plantdeepsea.ncpgr.cn/.

Bract suppression regulated by the miR156/529-SPLs-NL1-PLA1 module is required for the transition from vegetative to reproductive branching in rice

Reproductive transition of grasses is characterized by switching the pattern of lateral branches, featuring the suppression of outgrowth of the subtending leaves (bracts) and rapid formation of higher-order branches in the inflorescence (panicle). However, the molecular mechanisms underlying such changes remain largely unknown. Here, we show that bract suppression is required for the reproductive branching in rice. We identified a pathway involving the intrinsic time ruler microRNA156/529, their targets SQUAMOSA PROMOTER BINDING PROTEIN LIKE (SPL) genes, NECK LEAF1 (NL1), and PLASTOCHRON1 (PLA1), which regulates the bract outgrowth and thus affects the pattern switch between vegetative and reproductive branching. Suppression of the bract results in global reprogramming of transcriptome and chromatin accessibility following the reproductive transition, while these processes are largely dysregulated in the mutants of these genes. These discoveries contribute to our understanding of the dynamic plant architecture and provide novel insights for improving crop yields.

Discovering Cellular Mitochondrial Heteroplasmy Heterogeneity with Single Cell RNA and ATAC Sequencing

Next-generation sequencing technologies have revolutionised the study of biological systems by enabling the examination of a broad range of tissues. Its application to single-cell genomics has generated a dynamic and evolving field with a vast amount of research highlighting heterogeneity in transcriptional, genetic and epigenomic state between cells. However, compared to these aspects of cellular heterogeneity, relatively little has been gleaned from single-cell datasets regarding cellular mitochondrial heterogeneity. Single-cell sequencing techniques can provide coverage of the mitochondrial genome which allows researchers to probe heteroplasmies at the level of the single cell, and observe interactions with cellular function. In this review, we give an overview of two popular single-cell modalities-single-cell RNA sequencing and single-cell ATAC sequencing-whose throughput and widespread usage offers researchers the chance to probe heteroplasmy combined with cell state in detailed resolution across thousands of cells. After summarising these technologies in the context of mitochondrial research, we give an overview of recent methods which have used these approaches for discovering mitochondrial heterogeneity. We conclude by highlighting current limitations of these approaches and open problems for future consideration.

Fish-Ing for Enhancers in the Heart

Precise control of gene expression is crucial to ensure proper development and biological functioning of an organism. Enhancers are non-coding DNA elements which play an essential role in regulating gene expression. They contain specific sequence motifs serving as binding sites for transcription factors which interact with the basal transcription machinery at their target genes. Heart development is regulated by intricate gene regulatory network ensuring precise spatiotemporal gene expression program. Mutations affecting enhancers have been shown to result in devastating forms of congenital heart defect. Therefore, identifying enhancers implicated in heart biology and understanding their mechanism is key to improve diagnosis and therapeutic options. Despite their crucial role, enhancers are poorly studied, mainly due to a lack of reliable way to identify them and determine their function. Nevertheless, recent technological advances have allowed rapid progress in enhancer discovery. Model organisms such as the zebrafish have contributed significant insights into the genetics of heart development through enabling functional analyses of genes and their regulatory elements in vivo. Here, we summarize the current state of knowledge on heart enhancers gained through studies in model organisms, discuss various approaches to discover and study their function, and finally suggest methods that could further advance research in this field.

Integration of mRNA and miRNA Analysis Reveals the Molecular Mechanism of Cotton Response to Salt Stress

To identify the regulatory network of known and novel microRNAs (miRNAs) and their targets responding to salt stress, a combined analysis of mRNA libraries, small RNA libraries, and degradome libraries were performed. In this study, we used unique molecular identifiers (UMIs), which are more sensitive, accurate, and reproducible than traditional methods of sequencing, to quantify the number of molecules and correct for amplification bias. We identified a total of 312 cotton miRNAs using seedlings at 0, 1, 3, and 6 h after NaCl treatment, including 80 known ghr-miRNAs and 232 novel miRNAs and found 155 miRNAs that displayed significant differential expression under salt stress. Among them, fifty-nine differentially expressed miRNAs were simultaneously induced in two or three tissues, while 66, 11, and 19 were specifically expressed in the roots, leaves, and stems, respectively. It is indicated there were different populations of miRNAs against salt stress in roots, leaves and stems. 399 candidate targets of salt-induced miRNAs showed significant differential expression before and after salt treatment, and 72 targets of 25 miRNAs were verified by degradome sequencing data. Furthermore, the regulatory relationship of miRNA-target gene was validated experimentally via 5'RLM-RACE, proving our data reliability. Gene ontology and KEGG pathway analysis found that salt-responsive miRNA targets among the differentially expressed genes were significantly enriched, and mainly involved in response to the stimulus process and the plant hormone signal transduction pathway. Furthermore, the expression levels of newly identified miRNA mir1 and known miRNAs miR390 and miR393 gradually decreased when subjected to continuous salt stress, while overexpression of these miRNAs both increased sensitivity to salt stress. Those newly identified miRNAs and mRNA pairs were conducive to genetic engineering and better understanding the mechanisms responding to salt stress in cotton.