Method for predicting poplar root-specific cis-regulatory elements and modules at whole genome level

Abstract

The present application belongs to the technical field of bioinformatics, and particularly relates to a method for predicting poplar root-specific cis-regulatory elements and modules at the whole genome level. The present application discloses a method for predicting poplar root-specific cis-regulatory elements and modules at the whole genome level. For the first time, FIMO and tidymodels are combined to predict poplar tissue-specific cis-regulatory elements and modules at the whole genome level. Moreover, the promoter region is subjected to truncation analysis, and the TFBSs of the five promoter regions are respectively used to construct a random forest prediction model, so as to obtain TFBSs with high feature scores as candidate CREs. The candidate CREs are used to screen root-specific cis-regulatory modules, and experimental verification is performed. The method of the present application has the characteristics of rapidity and high throughput, and provides technical support for predicting plant root-specific cis-regulatory modules, and has practical value and wide application prospect.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of bioinformatics, and particularly relates to a method for predicting poplar root-specific cis-regulatory elements and modules at the whole genome level. BACKGROUND

[0002] A fundamental problem in biology is how to regulate the complex patterns of gene expression. At the core is the whole genome identification and characterization of cis-regulatory elements (CREs) and cis-regulatory modules (CRMs) that influence the expression of protein-coding and long non-coding RNA (lncRNA) genes. Here, a CRE is referred to as a single transcription factor binding site (TFBS), while a CRM is a collection of CREs, including promoters, transcriptional enhancers, silencers, and insulators. CRMs determine in which cell, at what time, and at what level a gene is expressed. In animals, chromatin accessibility and epigenomic mapping have made the identification of gene regulatory sequences routine. Precise positioning and characterization of CRMs and their target genes are challenging tasks; however, a variety of methods have proven effective.

[0003] One approach to identify unknown CRMs is to develop transgenic enhancer-trap lines, several lines with tissue-specific expression patterns have been identified in Arabidopsis, and enhancers have been successfully identified. While the above techniques are useful in dealing with small genome organisms that are easy to transform, for large genome organisms such as maize, quantitative trait locus (QTL) mapping and genetic fine mapping methods are used. The advent of high-throughput sequencing technology has greatly facilitated the identification of CRMs in the whole genome. In particular, epigenomic features, TF binding, and chromatin interactions have proven useful for CRM detection. Single-cell ATAC-seq is very effective for detecting CRMs at single-cell resolution. A major challenge at present is to find plant tissue- and environment-specific CRMs.

[0004] Machine learning (ML) methods are a promising alternative to inferential statistics that can be used to infer the cis-regulatory code from whole genome annotation data. In the human genome core promoters, the downstream core promoter region (DPR), DNA and RNA binding sequences, and other features have been determined by machine learning. In a recent plant study, a comprehensive study of transcription factor binding profiles predicted or inferred the cis-regulatory code behind transcriptional gene regulation. SUMMARY

[0005] The purpose of the present application is to use bioinformatics technology to screen poplar root tissue-specific CRMs from the whole genome in a high-throughput and rapid manner.

[0006] The present application provides the following technical solutions:

[0007] According to the poplar whole genome sequencing data, FIMO is used to detect poplar TFBSs elements in the whole genome, and on this basis, the distribution of TFBSs in the promoter of each gene is obtained, then the tidymodels is used to construct a random forest prediction model to obtain the feature score of TFBSs, and the TFBSs with high feature score are screened for predicting root-specific CRMs. The two bioinformatics tools are combined to complete the screening of poplar whole genome root-specific CREs. Because the arrangement and conservation of CREs are very complex, it is relatively difficult to accurately identify root-specific elements using most existing biological software. FIMO can search larger scale genomic data including whole genome sequence, and has higher accuracy compared with other software, and has the characteristics of high throughput; and the tidymodels is a machine learning integrated software, which can conveniently and efficiently construct and train the model.

[0008] The application comprises the following steps:

[0009] S101, using FIMO, taking the whole genome sequence as an input file, obtaining the distribution of poplar transcription factor binding sites (TFBSs) in the whole genome;

[0010] S102, using a self-programmed python script to count the distribution of TFBSs in the promoter of each gene;

[0011] S103, combining poplar RNA-seq data and tissue-specific parameters (τ), grouping the genes;

[0012] S104, using the tidymodels software, taking the distribution characteristics of TFBSs in the promoter of each gene as an input matrix, and the gene grouping in S103 as a label, to construct a random forest prediction model;

[0013] S105, adjusting the prediction model parameters to obtain high performance evaluation score and feature importance score;

[0014] S106, according to the feature importance score, screening the TFBSs with high feature score as candidate CREs;

[0015] S107, using the TFBSs with high feature score to screen the cis-regulatory modules, and performing experimental verification.

[0016] The application discloses a method for predicting poplar tissue-specific cis-regulatory elements and modules at the whole genome level by combining FIMO and tidymodels. The promoter region is subjected to truncation analysis, and TFBSs of five promoter regions are used to construct a random forest prediction model to obtain TFBSs with high feature scores as candidate CREs. The candidate CREs are used to screen root-specific cis-regulatory modules, and experimental verification is performed. BRIEF DESCRIPTION OF DRAWINGS

[0017] Figure 1 Fig. 1 is a schematic diagram of TFBSs distribution statistics of FIMO results of Populus trichocarpa;

[0018] Figure 2 Fig. 2 is a schematic diagram of learning input data of Populus trichocarpa obtained by using a python script;

[0019] Figure 3 Fig. 3 is a schematic diagram of TFBSs distribution statistics of FIMO results of Populus tomentosa;

[0020] Figure 4 Fig. 4 is a schematic diagram of learning input data of Populus tomentosa obtained by using a python script;

[0021] Figure 5 Fig. 5 is an expression pattern of PopRTS1 in 84K;

[0022] Figure 6 Fig. 6 is a TFBSs distribution in a PopRTS1 promoter;

[0023] Figure 7 Fig. 7 is a regulation of RTS on expression of iGUS gene in tobacco;

[0024] Figure 8 Fig. 8 is a flowchart of a method for predicting poplar root-specific cis-regulatory elements and modules at the whole genome level;

[0025] Figure 9 Fig. 9 is a python script for processing FIMO results. DETAILED DESCRIPTION

[0026] In order to further illustrate the application, the method for predicting poplar root-specific cis-regulatory elements and modules at the whole genome level provided by the application is described in detail below in combination with examples, but they should not be understood as limiting the protection scope of the application.

[0027] Example 1

[0028] Example 1: According to the process shown in Fig. 1, TFBSs in the whole genome of Populus trichocarpa are predicted and a model is constructed Figure 8

[0029] ​Selecting the Populus trichocarpa genome sequencing file, using FIMO software to identify the TFBSs elements at the whole genome level, and taking the results as the input file of machine learning, constructing a random forest model to predict the importance score of TFBSs.

[0030] Running steps:

[0031] 1. FIMO running

[0032] Although FIMO is composed of multiple perl language scripts, the inventors can run it with one command, taking the Populus trichocarpa genome as an example:

[0033] fimo [options] <motif file> <sequence file>

[0034] <motif file> is the PWM file of 379 Populus trichocarpa TFBSs (downloaded from PlantTFDB v5.0 database), with the file suffix meme (XXX.meme);

[0035] <sequence file> is the predicted DNA sequence in fasta format.

[0036] 2. FIMO output result processing

[0037] After running the output result, the gene start codon (ATG) upstream 3000bp and downstream 200bp are defined as the promoter, and the 3200bp promoter region is divided into +200 ~ -1000 (1200bp), +200 ~ -1500 (1700bp), +200 ~ 2000 (2200bp), +200 ~ -2500 (2700bp), +200 ~ -3000bp (3200bp) 5 regions; using python script ( Figure 9 ) to respectively count the distribution of each TFBS in each promoter in each region. The results (part) are shown in Figure 1 .

[0038] 3. FIMO output result processing

[0039] First, combined with RNA-seq data and tissue-specific parameters (τ), the genes are grouped, and the calculation method of τ is as follows:

[0040]

[0041] Where n is the number of tissues, x i is the average expression level of the gene in tissue i. is the expression profile component normalized by the maximum component value. The index varies from 0 to 1, where 0 means the expression of a gene is the same in all tissues, and 1 is tissue-specific. According to the annotated gene expression in the P. trichocarpa RNA-seq data, genes with no expression were determined as “silent genes”. The τ value was used to identify tissue-specific genes. Genes (τ = 1) were called “tissue-specific”, and genes (τ < 0.85, coefficient of variation CV < 0.30) were called “constitutively expressed”. The rest of the genes were labeled as “moderate frequency of expression”.

[0042] The python script (the script code of this step of the patent is written by the inventor) was used to preliminarily organize the promoter FIMO results of each region size, and to give each gene a corresponding classification label, to obtain the input data for machine learning, as shown in Figure 2

[0043] 4. Random forest model construction

[0044] The TFBSs distribution was used as a feature as input data in the ranger engine in tidymodels, 80% of which was used as training data and 20% as prediction data to construct a random forest model. At the same time, ten-fold validation was used to improve the performance of the model. The tune_grid function was used for hyperparameter selection. The performance indicators (F1-score, AUROC) of all 5 models were obtained (Table 1). Finally, the feature importance score of each input TFBSs in the trained model was extracted. The top 50 TFBSs with high feature scores were selected as candidate CREs for P. trichocarpa tissue specificity (Table 2).

[0045] Table 1 Performance evaluation indicators of the random forest model constructed in P. trichocarpa

[0046] Promoter size accuracy f-score-auc roc_auc f-score-auc 1200 bp 0.727 0.725167 0.592 0.721043 1700 bp 0.728 0.729523 0.606 0.728431 2200 bp 0.724 0.726733 0.583 0.724794 2700 bp 0.725 0.721085 0.58 0.719268 3200 bp 0.723 0.726294 0.578 0.725327

[0047] Table 2 Top 50 TFBSs with high feature scores in the 5 promoter regions of P. trichocarpa

[0048]

[0049]

[0050]

[0051]

[0052]

[0053] Example 2: Prediction of TFBSs in the whole genome of P. tomentosa and model construction

[0054] ​The Populus tomentosa genome sequencing file was selected, and FIMO software was used to identify the TFBSs element at the whole genome level of Populus tremuloides, and the results were used as the input file for machine learning to build a random forest model to predict the importance score of TFBSs.

[0055] Running steps:

[0056] 1. FIMO running

[0057] Although FIMO is composed of multiple perl language scripts, the inventors can run it with one command, for example, using the Populus tomentosa genome:

[0058] fimo [options] <motif file> <sequence file>

[0059] The <motif file> is a file of 379 PWMs of Populus tremuloides TFBSs (downloaded from the PlantTFDB v5.0 database), with the file suffix meme (XXX.meme);

[0060] The <sequence file> is the predicted DNA sequence in fasta format, which is the Populus tomentosa genome sequence here.

[0061] 2. FIMO output result processing

[0062] After running the output result, the start codon (ATG) upstream 3000bp and downstream 200bp of the gene were defined as the promoter, and the 3200bp promoter region was divided into +200~ -300, +200~ -500, +200~ -1000, +200~ -1500, +200~ 2000, +200~ -2500, +200~ -3000bp seven regions; the distribution of each TFBS in each promoter in each region was counted by using python script ( Figure 9 ) respectively. The results (part) are shown in Figure 3 .

[0063] 3. FIMO output result processing

[0064] First, the genes were grouped in combination with RNA-seq data and tissue-specific parameters (τ), and then according to the expression of the annotated genes in the Populus tomentosa RNA-seq data, the genes that were not expressed were determined as "silent genes". The τ value was used to identify tissue-specific genes. Genes (τ = 1) were called "tissue-specific", and genes (τ < 0.85, coefficient of variation CV < 0.30) were called "constitutive expression". The rest of the genes were labeled as "moderate frequency of expression".

[0065] The python script is used to preliminarily arrange each region size promoter FIMO result, and give each gene a corresponding classification label, and get the input data of machine learning, as shown in Figure 4

[0066] 4. Random forest model construction

[0067] Using the ranger engine in tidymodels, the TFBS distribution is used as a feature as input data, of which 80% is used as training data and 20% is used as prediction data to construct a random forest model. At the same time, ten-fold cross-validation is used to improve the performance of the model. Then the tune_grid function is used for hyperparameter selection. The performance indicators (F1-score, AUROC) of all 20 models are obtained (Table 3). Finally, the feature importance score of each input TFBS in the trained model is extracted. The top 50 TFBSs with high feature scores are selected as candidate CREs (Table 4).

[0068] Table 3 Performance evaluation index of random forest model constructed in Populus tomentosa

[0069]

[0070]

[0071] Table 4 Top 50 TFBSs with high feature scores in 5 promoter regions of Populus tomentosa

[0072]

[0073]

[0074]

[0075]

[0076]

[0077]

[0078] Example 3: Functional verification of TFBSs with high feature scores

[0079] Experimental steps:

[0080] 1. Select root-specifically expressed gene PopRTS1 to perform qRT-PCR to verify gene expression.

[0081] (1) Obtain the expression pattern of PopRTS1 gene

[0082] ​According to the transcriptome data of P. trichocarpa, the tissue specificity parameter (τ) of each gene was calculated, and the gene Potri.001G284300 with the largest τ value was selected. The sequence of the gene was obtained from the PopGenIE website (https: / / plantgenie.org / ), and the homologous gene PopRTS1 sequence was obtained in the silver gland poplar clone 84K by using blast software.

[0083] (2) Design qRT-PCR primers to explore the relative expression of PopRTS1

[0084] The qRT-PCR primers of PopRTS1 (qPopRSF / R) and Actin (internal reference gene) (qP18SF / R) (Table 5) were designed on the website (https: / / sg.idtdna.com / Scitools / Applications / RealTimePCR / ). The primers were synthesized by Sheng Wu Bioengineering (Shanghai) Co., Ltd. The mRNA of roots, stems and leaves was obtained by using the Tian Gen RNA prep Pure Plant Kit kit (DP432). The cDNA of the roots, stems and leaves of 84K tissue culture seedlings for one month was obtained by using the Nuowei reverse transcription kit (R101-02). According to the ABI StepOne Plus instrument manual, qRT-PCR was performed by using DNA binding dye SYBR green, and finally the relative expression abundance of PopRTS1 in roots, stems and leaves was calculated by using 2 –ΔΔCT The results proved that PopRTS1 was specifically expressed in roots as Potri.001G284300. Figure 5 ).

[0085] Table 5 All primers used in Example 3

[0086] Primers ID Sequence Usage SEQ ID No: qNActinR ACTTCCGGACATCTGAACCT Tobacco housekeeping gene 1 qNActinF GATCTTGCTGGTCGTGATCT Tobacco housekeeping gene 2 qP18SF GGCATGGAAGGTGATGCAGATC Poplar housekeeping gene 3 qP18SR CTGTGTCAAACAAGAACTTGTCC Poplar housekeeping gene 4 qGUSR CGTAAGGGTAATGCGAGGTAC Vector reporter gene 5 qGUSF TCTACTTTACTGGCTTTGGTCG Vector reporter gene 6 qHygF CTCCAGTCAATGACCGCTGT Vector resistance gene 7 qHygR CTCGGAGGGCGAAGAATCTC Vector resistance gene 8

[0087] 2, Predict the distribution of TFBSs in the PopRTS1 promoter.

[0088] The 1500bp upstream of the start codon (ATG) of PopRTS1 gene was selected as the promoter region, and the TFBSs in the promoter region were scanned by using FIMO software and 397 PWMs of P. trichocarpa TFBSs, and the results are as follows Figure 6 In A, 146 TFBSs were identified, which were enriched in the-1250bp-1500bp, -900bp-1050bp, -500bp-700bp and-50bp-400bp regions.

[0089] 3, Screening of CRMs by using high characteristic score TFBSs

[0090] First, the PWMs of 397 *Populus hairy fruit* TFBSs were narrowed down to 72 *Populus hairy fruit* tissue-specific candidate CREs, and a total of 57 TFBSs were identified, enriched in the regions of -100bp to -400bp and -500bp to -750bp. Figure 6 (B). Then, 98 tissue-specific candidate CREs of Populus tomentosa were identified in the promoter region, and a total of 51 TFBSs were identified, which were also enriched in the -100bp to -400bp and -500bp to -750bp regions. Figure 6 (C). Therefore, the cis-regulation module (CRMs) RTS ( Figure 6 (D).

[0091] 4. CRM functionality verification

[0092] (1) Acquisition of the 781bp promoter of the target fragment

[0093] First, the 84K genome was extracted using the Tiangen DNAsecure Plant Kit (DP320). Based on the distribution of promoter TFBSs, amplification primers RP (SEQ ID No: 9, 5'-CACAATCTTCCCTCTCCCCAT-3') and FP (SEQ ID No: 10, 5'-CTTTCAGTAGCCATAATCAAC-3') were designed and synthesized to obtain the -1 to -781bp promoter region. Full-length gold nanoparticles were then used for further analysis. The Cloning Kit ligates a 781bp promoter into a cloning vector for sequencing verification.

[0094] (2) Constructing expression carriers

[0095] First, the pCAMBIA1391 vector was modified to eliminate the promoter driving of the GUS gene (SEQ ID No: 11). A linear vector was obtained using SalI and BamHI restriction endonucleases. Next, ligation primers pFR (SEQ ID No: 12, 5'-aagcttggctgcaggCTTTCAGTAGCCATAATCAACAATC-3') and pFR (SEQ ID No: 13, 5'-aggatagtgggattgCTCGCTATTCTGGGATTTTTACACA-3') were designed using the website (https: / / soft.transgen.com.cn / index.php). A linear 781bp promoter containing a restriction endonuclease cantilever was obtained from the cloning vector using pFR and pFR. Finally, a transgold assay kit was used to... Seamless Cloning and Assembly Kit (CU201-02) will connect linear vector and linear promoter to get expression vector. The ligation product is transformed into E. coli competent cells, and the positive recombinant plasmid is screened by PCR. The recombinant plasmid DNA is extracted, and the Agrobacterium (GV3101) competent cells are transformed (liquid nitrogen freeze-thaw direct transfer method), and the recombinants are obtained by PCR detection.

[0096] (3) Cultivation of tobacco transgenic plants

[0097] ① Pre-culture: select 1-month-old tobacco tissue culture seedlings, take the tender leaves of tobacco aseptic seedlings with scissors at intervals of 0.5 cm, and keep the cut leaves from breaking to facilitate the next step of infection. The leaves are pre-cultured in the regeneration medium (MS+6-BA 1.0 mg / l+IBA 0.1 mg) for 1-2 days. When the cut part just starts to swell, Agrobacterium infection can be performed.

[0098] ② Agrobacterium culture: take the stored Agrobacterium from the -20°C refrigerator, and suck a small amount into YEB liquid medium or streak culture on a sterile operation table. Pick single colonies and inoculate into liquid YEB medium, add 50 mg / l Kan, 200 rpm, 28°C shake bacteria overnight. The next day, dilute the Agrobacterium culture in liquid MS medium; the second shaking needs 6-10 hours, and OD600 is 0.2-0.5 when it can be used for transformation.

[0099] ③ Infection: immerse the pre-cultured tobacco leaves in the Agrobacterium dilution liquid for 3 min on the clean bench, and then use sterile filter paper to absorb the excess bacterial liquid, and continue to place it in the original culture medium.

[0100] ④ Dark culture: dark culture at 25°C for 2-4 days, until white Agrobacterium is produced around the leaves.

[0101] ⑤ Selection culture: transfer the co-cultured leaves to the tobacco regeneration medium with the addition of 500 mg / L carbenicillin or 250 mg / L cephalosporin, and 50 mg / L kanamycin. Transfer to light under the conditions of 16 h day, 8 h night photoperiod and 25°C for selection culture.

[0102] ⑥ Subculture: after 2-3 weeks of selection culture, the differentiated resistant adventitious buds are transferred to the corresponding selection medium for subculture and expansion culture.

[0103] ⑦Rooting culture: the small seedlings which grow to more than 1cm were transferred to solid MS medium containing 300mg / L carbenicillin and 30mg / L kanamycin to induce rooting for about 2-3 weeks. The formation of roots can be observed, and then the regenerated whole plants were obtained. The rooting medium was 1 / 2MS medium containing 500mg / L cefotaxime and 100mg / L kanamycin.

[0104] (4) qRT-PCR identification of transgenic tobacco

[0105] ①Transgenic tobacco root, stem and leaf mRNA were obtained respectively by using Tian Gen RNA prep Pure Plant Kit kit (DP432).

[0106] ②Transgenic tobacco tissue culture seedling root, stem and leaf cDNA were obtained respectively by using Nuowei reverse transcription kit (R101-02).

[0107] ③According to the ABI StepOne Plus instrument manual, tobacco internal reference gene Actin was used as a control, qRT-PCR was performed by using DNA binding dye SYBR green, and finally 2 –ΔΔCT Method was used to calculate the relative non-expression of GUS and resistance gene in root, stem and leaf. The results proved that the -1-781bp promoter in PopRTS1 promoter sequence regulated the specific expression of GUS gene in root ( Figure 7 ).

[0108] Although the above embodiment has made a detailed description of the present application, it is only a part of the embodiments of the present application, but not all the embodiments. People can also obtain other embodiments according to the present embodiment without creativity, and these embodiments all belong to the protection scope of the present application.

Claims

1. A method for predicting tissue-specific cis-regulatory elements and modules in poplar root at the whole-genome level, characterized in that, Includes the following steps: S101. Using FIMO, with the whole genome sequence as input, obtain the distribution of transcription factor binding sites (TFBSs) in the poplar genome; step S101 includes: S1011. Obtain the PWM matrix of poplar TFBSs from the database website; S1012. Define the start codon ATG-3000bp ~ +200bp of each gene as the promoter region; S102. Statistically analyze the distribution pattern of TFBSs in the promoter of each gene using a self-written Python script; step S102 includes: S1021. Divide the 3200bp startup sub-region into 5 regions: upstream 200bp ~ downstream 1000bp, upstream 200bp ~ downstream 1500bp, upstream 200bp ~ downstream 2000bp, upstream 200bp ~ downstream 2500bp, and upstream 200bp ~ downstream 3000bp. S1022. Statistically analyze the distribution characteristics of TFBSs in the five promoter sub-regions respectively; S103. Combine RNA-seq data and tissue-specific parameter τ to group genes; S104. Using the tidymodels software, the distribution characteristics of TFBSs in the promoter of each gene are used as the input matrix, and the gene grouping in S103 is used as the label to construct a random forest prediction model. S105. Adjust the prediction model parameters to obtain the high-performance evaluation score and feature importance score; S106. Based on the feature importance scores, select TFBSs with high feature scores as candidate CREs; S107. Using TFBSs with high feature scores, cis-regulatory modules were screened and experimentally verified.

Images