# Comparative Transcriptomic Analysis of Root Cadmium Responses in Two Chinese Rice Cultivars Yuzhenxiang and Xiangwanxian 12

Authors: Shangdu Zhang; Zhenliang Luo; Xiang Wu; Bangzhi Shi; Huidan Jiang; Leliang Zhou; Lianyang Bai
Journal: Journal of Chemistry
 (2021)

---

## Abstract

Cadmium (Cd) pollution in paddy soil is an increasingly serious issue in rice production. It has been reported that there is a higher or lower grain Cd accumulation in the rice cultivars Yuzhenxiang (YZX) or Xiangwanxian 12 (XWX), respectively. To better manage the Cd pollution problem, the genes that might play vital roles in governing the difference in root Cd responses between these two rice cultivars were examined. In this study, the results of RNA sequencing (RNA-seq) showed that there were 341 and 161 differentially expressed genes in the roots of YZX and XWX after Cd exposure, respectively. Among these genes, 7 genes, such as Os06g0196300 (OsJ_019618), Os07g0570700 (OsJ_24808), ADI1, GDCSH, HSFB2C, PEX11-4, and CLPB1, possessed higher degree nodes with each other, through interaction analysis by the STRING (search tool for the retrieval of interacting genes/proteins) software, suggesting that they might play vital roles in Cd response. Based on GO enrichment analysis, 41 differently expressed genes after Cd treatment in YZX or XWX were identified to be related to Cd response. Through comparative transcriptomic analysis, 257 genes might be associated with the root Cd response difference between YZX and XWX. Furthermore, we supposed that ADI1, CFBP1, PEX11-4, OsJ_019618, OsJ_24808, GDCSH, CLPB1, LAC6, and WNK3 might be implicated in Cd response based on the combined analysis of RT-qPCR, interaction, and GO annotation analysis. In conclusion, the numerous genes that might be related to Cd stress response and root Cd response difference between YZX and XWX at the booting stage may be of benefit for the development of rice varieties with low Cd consumption.


---

## Body

## 1. Introduction

Rice (Oryza sativa L.) is the primary source of calorie intake in plenty of countries, including China, and the staple food of over half of the world’s population [1, 2]. To meet the food demand of increasing world populations, it is imperative to improve the rice yield [3]. However, the sustainability and productivity of the rice production system are constantly threatened by multiple factors, such as environmental pollution, climate change, and excessive use of pesticides and fertilizers [2].Heavy metal pollution is the predominant type in all soil contaminants, and cadmium (Cd) contamination ranks the first among heavy metal contaminants in China [4, 5]. Cd contamination, mainly caused by mining, e-waste, and overuse of nitrogen and phosphate fertilizers, has been found in large-scale agricultural soil in China [6, 7]. Cd is readily absorbed by rice plants and easily transferred to food chains [8]. Moreover, it is often accumulated in rice grains and human bodies due to its long half-life of up to 25–30 years [8, 9]. Excessive Cd exposure not only can inhibit the growth of rice plants but also can decrease the quality/nutrients/yields of rice grains [7]. Additionally, the intake of Cd-contaminated rice can markedly increase the risks of multiple diseases, such as cancer, osteoporosis, and liver/kidney injury, as well as nervous system diseases [9–11]. To reduce the potential harm of Cd on human health, it is imperative to screen out the genes that play vital roles in Cd response in rice plants [7].It has been reported that there is a notable difference in Cd content in the grains of different rice cultivars, which might be caused by the genotypic and environmental diversities among rice cultivars [12, 13]. More specifically, rice Cd concentration is closely associated with heading time, soil pH, mutation/dysregulation/diversity of multiple genes, or quantitative trait locus related to Cd response (absorption, translocation, and accumulation) among different rice subspecies and cultivars [12, 13]. Over the past decades, high-throughput RNA sequencing (RNA-seq) in combination with function annotation/enrichment and bioinformatics prediction analysis has been widely used to decipher or speculate the essential genes/biological processes/pathways under different conditions at the transcriptomic level and gain a deep and comprehensive understanding of molecular basis underlying the phenotypic/biological differences in plants including rice [14–16].In this study, RNA-seq-based transcriptomic analysis of root samples of YZX and XWX with or without Cd treatment at the booting stage was performed to identify key genes in root Cd response in Yuzhenxiang (YZX, a high Cd accumulation rice cultivar) and Xiangwanxian 12 (XWX, a low Cd accumulation rice cultivar) and genes associated with root Cd accumulation difference between YZX and XWX. RNA-seq results showed that 341 and 161 transcripts were differentially expressed in the roots of YZX or XWX rice after Cd treatment. Moreover, the protein-protein interaction (PPI) networks of dysregulated transcripts in response to Cd exposure in YZX or XWX were established by the STRING database and 7 genes that might play vital roles in the response to Cd pollution were screened out based on the node degrees of proteins in the PPI networks. Additionally, 41 transcripts that function as crucial players in Cd response were filtered out based on GO annotation analysis. Furthermore, 257 transcripts that might be associated with the difference in root Cd response between YZX or XWX were screened out by comparative transcriptomic analysis. Also, we further examined the expression patterns of 10 genes (ferredoxin-1 (ADI1), fructose-1, 6-bisphosphatase (CFBP1), glycine cleavage system H protein (GDCSH), laccase-6 (LAC6), peroxiredoxin Q (Os06g0196300), ribosome-recycling factor (Os07g0570700), peroxisomal membrane protein 11–4 (PEX11-4), probable serine/threonine-protein kinase WNK3, chaperone protein CLPB1, and heat stress transcription factor B-2c (HSFB2C)) by RT-qPCR assay and filtered some key genes related to Cd stimulation based on RNA-seq/RT-qPCR/interaction/GO enrichment/GO annotation analysis together with previous study outcomes.


## 2. Materials and Methods

### 2.1. Plant Materials and Treatment

TwoIndica rice cultivars (i.e., YZX and XWX) were cultivated in the LT-36VL climatic growth chamber (Percival, Perry, IA, USA) at the Hunan Agricultural Biotechnology Research Institute (Changsha, China). In the preliminary tests, YZX and XWX rice at 4 different growth phases (i.e., seedling, tillering, booting, and grain-filling stages) were cultivated in a hydroponic system containing 2 mg/L of Cd under low-temperature (15–20°C), moderate-temperature (22–27°C), and high-temperature (30–35°C) condition for 48 h. Cd concentration was determined according to the instructions in the Chinese National Standard GB 5009.15–2014.


### 2.2. RNA Sequencing

In the RNA-seq and real-time quantitative reverse transcription PCR (RT-qPCR) validation experiments, YZX and XWX rice were grown at 22–27°C (moderate temperature) in the hydroponic system with or without 2 mg/L of Cd exposure for 48 h. Next, rice root samples were collected. The samples were divided into 4 groups: AMD group (YZX group with Cd treatment), AMK group (YZX group without Cd treatment), BMD group (XWX group with Cd addition), and BMK group (XWX group without Cd exposure). Each group contains 3 root samples.RNA was extracted from the root samples using the TRIzol reagent (Thermo Scientific, Waltham, MA, USA) according to the protocols of the manufacturer. The concentration, purity, and quality of RNA were examined by NanoDrop 2000 spectrophotometer (Thermo Scientific). RNA sequencing was carried out as previously described [17]. The cDNA library was constructed using the TruSeq RNA sample preparation kit (Illumina, San Diego, CA, USA) following the protocols of the manufacturer. The mRNA with poly-A structure was enriched by Poly-T oligo-attached magnetic beads (Illumina). The cDNA library was sequenced using the PE150 strategy by Biomarker Technologies Co., Ltd. (Beijing, China). Raw sequencing data were processed and filtered into high-quality clean data. Next, the clean data were aligned to the reference genome of Oryza sativa L. subsp. Indica using HISAT2 and then assembled and quantified using StringTie. The genome version used in the sequence assembly was ASM465v1. (http://plants.ensembl.org/Oryza_indica/Info/Index).


### 2.3. Differential Expression Analysis

Differential expression analysis was performed using the DESeq R package. Genes were defined to be differentially expressed at the fold change ≥2 and FDR <0.01.


### 2.4. Gene Enrichment Analysis

KEGG enrichment analysis was performed using the KOBAS3 database (http://kobas.cbi.pku.edu.cn/kobas3/genelist/). P<0.05 represented that the terms were significantly enriched.


### 2.5. Interaction Network Construction

The interaction networks of gene-coding proteins were constructed by the STRING database (https://www.string-db.org/).


### 2.6. RT-qPCR Assay

RNA was reversely transcribed into first-strand cDNA without genomic DNA contamination using TransScript One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen Biotech Co., Ltd, Beijing, China). Briefly, RNA was coincubated with Oligo (dT) 18 primer for 5 min at 65°C and then immediately placed on an ice bath for 2 min. Then, the reaction mixture was coincubated with TS Reaction Mix, TransScript RT/RI Enzyme Mix, and gDNA Remover at 42°C for 30 min and then terminated at 85°C for 5 min. Next, cDNA was amplified and quantified using the ChamQ Universal SYBR qPCR Master Mix (Vazyme Biotech Co., Ltd., Nanjing, China) and specific quantitative primers under the following thermocycling conditions: 95°C for 3 min, 40 cycles of 95°C for 10 s, and 60°C for 30 s. 18sRNA served as the housekeeping gene to normalize the expression of other genes. RT-qPCR was performed on Roche LightCycler 480 II Real-Time PCR Detection System (Basel, Sweden). The quantitative primer sequences are presented in Table1. The RNA samples in the RT-qPCR assay were identical to those used in RNA-seq experiments.Table 1
Quantitative primer sequences of 10 interested genes.
Gene IDPrimer namePrimer (5′-3′)18s rRNA18s rRNA-F1CTACGTCCCTGCCCTTTCTACA18s rRNA-R1ACACTTCACCGGACCATTCAABGIOSGA029403HSFB2C–F1ACAACTTCTCCAGCTTCGTGHSFB2C-R1ACCTTCCGGCGGTGTATATCBGIOSGA016829PEX11-4-F1TCGTGAGCTGTTCCAACTCTPEX11-4-R1TGCCAGTGAGCTGTTCAGTABGIOSGA017714CLPB1-F1GTTCGAAGAACGGCTCAAGGCLPB1-R1ACGAGGTGTATCTCGTCGATBGIOSGA015634WNK3–F1CGTCGACTTCATGTCATCGGWNK3-R1TGCAAATGCTGCTCCTCTTCBGIOSGA033278GDCSH-F1AGGACGGGTGGATGATCAAGGDCSH-R1AGAGCCTAGTGAGCGTCTTCBGIOSGA004865CFBP1–F1GGGAAGTATTGCGTGTGCTTCFBP1-R1AGGCTGTAACACGTCCTCAABGIOSGA027798ADI1-F1TCATCGAGACCCACAAGGAGADI1-R1AGCATAGGGACGACGACATTBGIOSGA004807LAC6-F1TGTCAGCCAGGGAATACGAGLAC6-R1GGACTGTCACGTTGTATGGCBGIOSGA025975Os07g0570700-F1GTTCCAGAACGCCCAAACATOs07g0570700-R1GGCATGCCTCAATACTGCTCBGIOSGA021734Os06g0196300-F1ACAAGGTGAGGAAGGAGTGGOs06g0196300-R1GGAGGATCTTGAGGGTCTCG


### 2.7. Statistical Analysis

Data analysis was performed using GraphPad Prism 7 software (GraphPad Software, Inc., San Diego, CA, USA). Results were displayed as mean ± standard deviation. The difference among groups was examined using one-way ANOVA and Dunnett post hoc test. Statistically significant differences were set atP<0.05.


### 2.8. Accession Numbers

The raw data of RNA-seq reads were deposited in the National Center for Biotechnology Information (NCBI) database under the accession number (Submission ID: SUB10609471, BioProject ID: PRJNA777353). Biosample accessions were as follows: SAMN22851327, SAMN22851328, SAMN22851329, SAMN22851330, SAMN22851331, SAMN22851332, SAMN22851333, SAMN22851334, SAMN22851335, SAMN22851336, SAMN22851337, and SAMN22851338.


## 2.1. Plant Materials and Treatment

TwoIndica rice cultivars (i.e., YZX and XWX) were cultivated in the LT-36VL climatic growth chamber (Percival, Perry, IA, USA) at the Hunan Agricultural Biotechnology Research Institute (Changsha, China). In the preliminary tests, YZX and XWX rice at 4 different growth phases (i.e., seedling, tillering, booting, and grain-filling stages) were cultivated in a hydroponic system containing 2 mg/L of Cd under low-temperature (15–20°C), moderate-temperature (22–27°C), and high-temperature (30–35°C) condition for 48 h. Cd concentration was determined according to the instructions in the Chinese National Standard GB 5009.15–2014.


## 2.2. RNA Sequencing

In the RNA-seq and real-time quantitative reverse transcription PCR (RT-qPCR) validation experiments, YZX and XWX rice were grown at 22–27°C (moderate temperature) in the hydroponic system with or without 2 mg/L of Cd exposure for 48 h. Next, rice root samples were collected. The samples were divided into 4 groups: AMD group (YZX group with Cd treatment), AMK group (YZX group without Cd treatment), BMD group (XWX group with Cd addition), and BMK group (XWX group without Cd exposure). Each group contains 3 root samples.RNA was extracted from the root samples using the TRIzol reagent (Thermo Scientific, Waltham, MA, USA) according to the protocols of the manufacturer. The concentration, purity, and quality of RNA were examined by NanoDrop 2000 spectrophotometer (Thermo Scientific). RNA sequencing was carried out as previously described [17]. The cDNA library was constructed using the TruSeq RNA sample preparation kit (Illumina, San Diego, CA, USA) following the protocols of the manufacturer. The mRNA with poly-A structure was enriched by Poly-T oligo-attached magnetic beads (Illumina). The cDNA library was sequenced using the PE150 strategy by Biomarker Technologies Co., Ltd. (Beijing, China). Raw sequencing data were processed and filtered into high-quality clean data. Next, the clean data were aligned to the reference genome of Oryza sativa L. subsp. Indica using HISAT2 and then assembled and quantified using StringTie. The genome version used in the sequence assembly was ASM465v1. (http://plants.ensembl.org/Oryza_indica/Info/Index).


## 2.3. Differential Expression Analysis

Differential expression analysis was performed using the DESeq R package. Genes were defined to be differentially expressed at the fold change ≥2 and FDR <0.01.


## 2.4. Gene Enrichment Analysis

KEGG enrichment analysis was performed using the KOBAS3 database (http://kobas.cbi.pku.edu.cn/kobas3/genelist/). P<0.05 represented that the terms were significantly enriched.


## 2.5. Interaction Network Construction

The interaction networks of gene-coding proteins were constructed by the STRING database (https://www.string-db.org/).


## 2.6. RT-qPCR Assay

RNA was reversely transcribed into first-strand cDNA without genomic DNA contamination using TransScript One-Step gDNA Removal and cDNA Synthesis SuperMix (TransGen Biotech Co., Ltd, Beijing, China). Briefly, RNA was coincubated with Oligo (dT) 18 primer for 5 min at 65°C and then immediately placed on an ice bath for 2 min. Then, the reaction mixture was coincubated with TS Reaction Mix, TransScript RT/RI Enzyme Mix, and gDNA Remover at 42°C for 30 min and then terminated at 85°C for 5 min. Next, cDNA was amplified and quantified using the ChamQ Universal SYBR qPCR Master Mix (Vazyme Biotech Co., Ltd., Nanjing, China) and specific quantitative primers under the following thermocycling conditions: 95°C for 3 min, 40 cycles of 95°C for 10 s, and 60°C for 30 s. 18sRNA served as the housekeeping gene to normalize the expression of other genes. RT-qPCR was performed on Roche LightCycler 480 II Real-Time PCR Detection System (Basel, Sweden). The quantitative primer sequences are presented in Table1. The RNA samples in the RT-qPCR assay were identical to those used in RNA-seq experiments.Table 1
Quantitative primer sequences of 10 interested genes.
Gene IDPrimer namePrimer (5′-3′)18s rRNA18s rRNA-F1CTACGTCCCTGCCCTTTCTACA18s rRNA-R1ACACTTCACCGGACCATTCAABGIOSGA029403HSFB2C–F1ACAACTTCTCCAGCTTCGTGHSFB2C-R1ACCTTCCGGCGGTGTATATCBGIOSGA016829PEX11-4-F1TCGTGAGCTGTTCCAACTCTPEX11-4-R1TGCCAGTGAGCTGTTCAGTABGIOSGA017714CLPB1-F1GTTCGAAGAACGGCTCAAGGCLPB1-R1ACGAGGTGTATCTCGTCGATBGIOSGA015634WNK3–F1CGTCGACTTCATGTCATCGGWNK3-R1TGCAAATGCTGCTCCTCTTCBGIOSGA033278GDCSH-F1AGGACGGGTGGATGATCAAGGDCSH-R1AGAGCCTAGTGAGCGTCTTCBGIOSGA004865CFBP1–F1GGGAAGTATTGCGTGTGCTTCFBP1-R1AGGCTGTAACACGTCCTCAABGIOSGA027798ADI1-F1TCATCGAGACCCACAAGGAGADI1-R1AGCATAGGGACGACGACATTBGIOSGA004807LAC6-F1TGTCAGCCAGGGAATACGAGLAC6-R1GGACTGTCACGTTGTATGGCBGIOSGA025975Os07g0570700-F1GTTCCAGAACGCCCAAACATOs07g0570700-R1GGCATGCCTCAATACTGCTCBGIOSGA021734Os06g0196300-F1ACAAGGTGAGGAAGGAGTGGOs06g0196300-R1GGAGGATCTTGAGGGTCTCG


## 2.7. Statistical Analysis

Data analysis was performed using GraphPad Prism 7 software (GraphPad Software, Inc., San Diego, CA, USA). Results were displayed as mean ± standard deviation. The difference among groups was examined using one-way ANOVA and Dunnett post hoc test. Statistically significant differences were set atP<0.05.


## 2.8. Accession Numbers

The raw data of RNA-seq reads were deposited in the National Center for Biotechnology Information (NCBI) database under the accession number (Submission ID: SUB10609471, BioProject ID: PRJNA777353). Biosample accessions were as follows: SAMN22851327, SAMN22851328, SAMN22851329, SAMN22851330, SAMN22851331, SAMN22851332, SAMN22851333, SAMN22851334, SAMN22851335, SAMN22851336, SAMN22851337, and SAMN22851338.


## 3. Results

### 3.1. Cd Accumulation Difference in the Roots of YZX and XWX at Four Different Growth Stages

Previous studies have reported that there is a higher Cd accumulation in the first node under panicle, panicle node, and grains of YZX cultivar compared to XWX cultivar [16, 18]. We supposed that the difference in grain Cd accumulation in these two rice cultivars might be caused by the discrepancy of root Cd absorptive ability. Moreover, root Cd uptake capacity has been found to be varied at different growth temperatures and growth stages in rice [19–21]. To identify the optimal conditions under which the difference of Cd concentration in the roots of YZX and XWX was the most prominent, the concentrations of Cd in the roots of these two rice cultivars at the seedling/tillering/booting/grain-filling stage under low/moderate/high temperature and 2 mg/L of Cd stress were measured. Results showed that there was a noticeable difference in root Cd content at the booting stage when these two rice cultivars were cultured under the moderate-temperature condition (Figure 1). Hence, YZX and XWX rice were cultivated at moderate temperature and 2 mg/L of Cd stress for 48 h to explore the gene expression alterations in response to Cd treatment in the following experiments.Figure 1
Cd accumulation difference in the roots of YZX and XWX at 4 different growth stages under different temperatures and the same Cd stress. YZX and XWX rice at different growth stages (i.e., seedling, tillering, booting, and grain-filling stages) were cultured for 48 h under different temperatures (low temperature (15–20°C), moderate temperature (22–27°C), or high temperature (30–35°C)) and 2 mg/L of Cd stress. Next, Cd concentrations in the roots of these two rice cultivars were determined.


### 3.2. Identification of Differentially Expressed Transcripts after Cd Treatment in YZX and XWX

In this study, RNA-seq technology was used to identify Cd response-related transcripts in YZX and XWX and investigate genes that might play vital roles in governing the difference in root Cd accumulation between YZX and XWX. The statistics information of clean data after filtration of RNA-seq raw data is shown in Table2. The GC content of clean data was approximately 50% (Table 2). Moreover, more than 94% of bases have a recognition accuracy of over 99.9% in the clean data (Table 2). Differential expression analysis showed that 341 transcripts were differentially expressed (206 downregulated and 135 upregulated) in the AMD vs. AMK group (Supplementary Figure 1(a) and Supplementary Table 1). Also, 75 downregulated transcripts and 86 upregulated transcripts were identified in the BMD vs. BMK group (Supplementary Figure 1(b) and Supplementary Table 2). Moreover, 18 dysregulated transcripts (i.e., BGIOSGA022091, BGIOSGA007908, BGIOSGA007302, BGIOSGA029201, BGIOSGA022060, BGIOSGA023662, BGIOSGA026675, BGIOSGA001422, BGIOSGA025476, BGIOSGA030712, BGIOSGA004724, BGIOSGA018872, BGIOSGA010590, BGIOSGA001025, BGIOSGA016478, BGIOSGA020263, BGIOSGA001117, and BGIOSGA000917) in the AMD vs. AMK group and 6 dysregulated transcripts (i.e., BGIOSGA004809, BGIOSGA030712, BGIOSGA034657, BGIOSGA012380, BGIOSGA018872, and BGIOSGA020584) in the BMD vs. BMK group were identified to be involved in inorganic ion transport and metabolism based on COG_class and KOG_class annotation analyses in Supplementary Tables 1 and 2. These transcripts might play vital roles in governing the difference in root Cd accumulation in YZX and XWX.Table 2
The statistics information of clean data.
Sample no.Clean readsClean basesGC content (%)% ≥Q30AMD-128,320,9218,469,123,35452.5194.89AMD-224,166,3847,224,223,02052.6295.04AMD-323,919,6997,147,482,54452.1695.05AMK-125,675,7837,673,764,99452.9994.99AMK-224,904,1077,449,500,70651.7595.23AMK-323,973,6977,172,841,03252.7994.77BMD-130,635,0409,159,932,34053.2896.41BMD-221,090,0586,311,270,66851.6494.84BMD-323,708,8487,088,906,52251.6794.81BMK-125,223,0397,541,950,04452.2894.88BMK-237,274,48111,138,840,83253.1396.36BMK-327,498,8598,216,506,79852.4996.44Notes. Clean reads: the total number of paired-end reads in clean data; clean bases: the total number of bases in clean data; GC content: the percentages of G and C bases in clean data; ≥Q30%: the percentage of bases with the quality value ≥ 30 (base identification accuracy ≥ 99.9%).


### 3.3. Construction of PPI Networks of Dysregulated Genes in the AMD vs. AMK and BMD vs. BMK Groups

Next, the STRING database was used to investigate the relationships of proteins encoded by the abovementioned differentially expressed genes. The interaction networks of the dysregulated genes in the AMD vs. AMK group are presented in Supplementary Figure2(a) and Supplementary Table 3. Similarly, the interaction relationships of the differentially expressed genes in the BMD vs. BMK group are shown in Supplementary Figure 2(b) and Supplementary Table 4. Additionally, the node degree (number of interacted proteins) of each protein in the PPI network of Supplementary Figures 2(a) and 2(b) was analyzed and is displayed in Supplementary Tables 5 and 6, respectively. Genes with greater node degrees might play central roles in the response to Cd. Among genes with higher node degrees, 7 genes (Os06g0196300 (OsJ_019618), Os07g0570700 (OsJ_24808), ADI1, GDCSH, HSFB2C, PEX11-4, and CLPB1) were screened out for further RT-qPCR validation. Moreover, interaction analysis of dysregulated genes in the AMD vs. AMK group revealed that there was a complex interaction among OsJ_019618, OsJ_24808, ADI1, GDCSH, and PEX11-4. For instance, OsJ_019618 could interact with OsJ_24808, ADI1, GDCSH, or PEX11-4 (Supplementary Tables 3 and 5). ADI1 could interact with OsJ_019618, GDCSH, or PEX11-4 (Supplementary Tables 3 and 5). Moreover, there was a potential interaction between OsJ_24808 and OsJ_019618 or GDCSH (Supplementary Tables 3 and 5).


### 3.4. Comparative Analysis of Cd-Responsive Genes in XWX and YZX

GO annotation analysis of the 341 differentially expressed transcripts in the AMD vs. AMK group showed that 37 transcripts were involved in the response to Cd (Supplementary Table7). Moreover, 23 of the 161 differentially expressed transcripts in the BMD vs. BMK group were identified to be implicated in Cd response based on GO annotation analysis of dysregulated transcripts in the BMD vs. BMK group (Supplementary Table 8). Among these Cd-responsive transcripts, 19 common elements were found in AMD vs. AMK and BMD vs. BMK groups (Figure 2(a) and Supplementary Table 9). These common transcripts had the same alteration trends in the AMD vs. AMK and BMD vs. BMK groups (Supplementary Table 9). And 18 or 4 transcripts were found to be differentially expressed only in the AMD vs. AMK group or BMD vs. BMK group, respectively (Figure 2(a) and Supplementary Table 9). These 22 transcripts might have crucial roles in regulating the difference in root Cd accumulation between YZX and XWX at the booting stage. Additionally, we noticed that 43 transcripts were differentially expressed only in the BMD vs. BMK group, but not notably altered in the AMD vs. AMK group (Supplementary Table 10). Also, 214 transcripts were found to be markedly upregulated or downregulated in the AMD vs. AMK group, whereas the difference in the expression of these 214 transcripts was not significant in the BMD vs. BMK group (Supplementary Table 11). The 257 transcripts in Supplementary Tables 10 and 11 are integrated into Supplementary Table 12. We believed that these transcripts might also be related to the difference in root Cd accumulation between XWX and YZX. Also, transcripts with notably different upregulated or downregulated ratios in the BMD vs. BMK and AMD vs. AMK groups might be implicated in root Cd accumulation difference between two rice cultivars, which were not analyzed in our project. KEGG enrichment analysis showed that the transcripts in Supplementary Table 12 mainly participated in the regulation of metabolic pathways, biosynthesis of secondary metabolites, carbon fixation in photosynthetic organisms, carbon metabolism, glutathione metabolism, fructose and mannose metabolism, glycolysis/gluconeogenesis, starch and sucrose metabolism, and glycerophospholipid metabolism (Supplementary Table 13 and Figure 2(b)). For instance, 23, 16, 5, 7, or 5 differentially expressed genes were significantly enriched in KEGG pathways related to metabolism, biosynthesis of secondary metabolites, carbon fixation in photosynthetic organisms, carbon metabolism, or glutathione metabolism, respectively (Supplementary Table 13). The top 20 KEGG pathways are shown in Figure 2(b).Figure 2
Comparative analysis of Cd-responsive genes in XWX and YZX. (a) Venn analysis for transcripts in Supplementary Tables7 and 8. AMD vs. AMK list: differentially expressed transcripts related to Cd response in the AMD vs. AMK group based on GO annotation analysis of dysregulated transcripts in the AMD vs. AMK group. BMD vs. BMK list: differentially expressed transcripts related to Cd response in the BMD vs. BMK group based on GO annotation analysis of dysregulated transcripts in the BMD vs. BMK group. Upper subfigure: Venn diagrams showing a number of differentially expressed transcripts related to Cd response in the AMD vs. AMK and BMD vs. BMK groups, as well as the number of overlapping genes. Middle subfigure: number of transcripts in the AMD vs. AMK and BMD vs. BMK lists. Bottom subfigure: number of specific (1) and shared (2) transcripts in the AMD vs. AMK and BMD vs. BMK lists. (b) The top 20 KEGG pathways enriched by the transcripts in Supplementary Table 12. The circle size denotes the number of transcripts in corresponding KEGG pathway terms. The circle color represents the “−log10 (P value)”.
(a)(b)


### 3.5. Expression Analysis of 10 Interested Genes by RT-qPCR Assay

Next, 10 genes including the 7 abovementioned genes with greater node degrees in the PPI networks, 1 interested gene (LAC6), and 2 genes related to heavy metal stress tolerance (WNK3 and CFBP1) (detailed information is shown in the Discussion section) were selected for further RT-qPCR validation. RT-qPCR results showed that the expression levels of ADI1, CFBP1, GDCSH, LAC6, Os06g0196300, Os07g0570700, PEX11-4, and WNK3 were notably downregulated in both AMD vs. AMK and BMD vs. BMK groups (Figure3(a), 3(b), 3(d), and 3(f)–3(j)). CLPB1 and HSFB2C were highly expressed in the AMD vs. AMK group, but not notably changed in the BMD vs. BMK group (Figures 3(c) and 3(e)). Comparative analysis of RT-qPCR and RNA-seq outcomes disclosed that the alteration trends of these 10 genes were consistent in response to Cd treatment in YZX and XWX rice cultivars, although some differences were not statistically significant (Figures 3(a)–3(j)).Figure 3
Comparative analysis of RT-qPCR and RNA-seq outcomes of expression patterns of 10 interested genes (i.e., PEX11-4 (a), WNK3 (b), CLPB1 (c), GDCSH (d), HSFB2C (e), LAC6 (f), Os06g0196300 (g), Os07g0570700 (h), ADI1 (i), and CFBP1 (j)) in roots of XWX and YZX.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)


### 3.6. Interaction and Enrichment Analysis of These 10 Genes

The relationships of these 10 genes were further deciphered by PPI analysis via the STRING database. Results suggested that 7 gene-encoded proteins (ADI1, CFBP1, PEX11-4, OsJ_019618, OsJ_24808, GDCSH, and CLPB1) could form a potential regulatory network (Figure4). The interaction analysis of these 10 proteins is shown in sheet 1 of Supplementary Table 14. Moreover, GO enrichment analysis revealed that 7 genes (LAC6, CFBP1, OsJ_019618, WNK3, OsJ_24808, ADI1, and GDCSH) were enriched in the metabolic process (Supplementary Table 14, sheet 2). Additionally, 6 genes (LAC6, CFBP1, CLPB1, WNK3, OsJ_24808, and ADI1) had a potential ion binding activity (Supplementary Table 14, sheet 3). Furthermore, annotation analysis revealed that 3 genes (LAC6, OsJ_019618, and PEX11-4) were involved in the regulation of the response to oxidative stress (Supplementary Table 14, sheet 4).Figure 4
Interaction network of 10 interested genes. The yellow line represents that the potential interaction between proteins was obtained by text mining. The black line denotes that the potential interaction between proteins was identified by coexpression analysis. The green line means that the potential interaction between proteins was predicted based on gene neighborhood relationships. PPI analysis suggested that PEX11-4 could interact with CFBP1 (FBP1), ADI1, and Os06g0196300 (OsJ_019618). There was a potential interaction between CFBP1 (FBP1) and PEX11-4/ADI1/Os06g0196300 (OsJ_019618)/GDCSH. GDCSH coud interact with CFBP1 (FBP1)/ADI1/Os06g0196300 (OsJ_019618)/Os07g0570700 (OsJ_24808). Os07g0570700 (OsJ_24808) could interact with GDCSH/Os06g0196300 (OsJ_019618). Os06g0196300 (OsJ_019618) could interact with ADI1/PEX11-4/CFBP1 (FBP1)/GDCSH/Os07g0570700 (OsJ_24808)/CLPB1. ADI1 could interact with PEX11-4/CFBP1 (FBP1)/GDCSH/Os06g0196300 (OsJ_019618). CLPB1 could interact with Os06g0196300 (OsJ_019618).


### 3.7. Correlation Analysis of These 10 Genes and Cd Response

As shown in Supplementary Table5, OsJ_019618, OsJ_24808, ADI1, GDCSH, PEX11-4, CLPB1, and HSFB2C could interact with 59, 44, 36, 29, 17, 6, and 6 genes in the interaction network of dysregulated genes in the AMD vs. AMK group, respectively. And CLPB1, PEX11-4, and HSFB2C could interact with 6, 5, and 3 genes in the BMD vs. BMK group, respectively (Supplementary Table 6). To explore the association of these genes and Cd response, Venn analysis of genes that could interact with the 10 abovementioned genes and GO-annotated Cd-related genes in Supplementary Tables 7 and 8 was performed. Results showed that OsJ_019618 could interact with 4 genes related to the Cd response (i.e., APX2, HCF136, GLN2, and ALDP). ADI1 could interact with 2 Cd-related genes (i.e., ALDP and GLN2). OsJ_24808, GDCSH, or CLPB1 could interact with Cd-responsive gene TPI, GLN2, or HSP17.4, respectively. These data further suggested the close association of these genes and Cd response.


## 3.1. Cd Accumulation Difference in the Roots of YZX and XWX at Four Different Growth Stages

Previous studies have reported that there is a higher Cd accumulation in the first node under panicle, panicle node, and grains of YZX cultivar compared to XWX cultivar [16, 18]. We supposed that the difference in grain Cd accumulation in these two rice cultivars might be caused by the discrepancy of root Cd absorptive ability. Moreover, root Cd uptake capacity has been found to be varied at different growth temperatures and growth stages in rice [19–21]. To identify the optimal conditions under which the difference of Cd concentration in the roots of YZX and XWX was the most prominent, the concentrations of Cd in the roots of these two rice cultivars at the seedling/tillering/booting/grain-filling stage under low/moderate/high temperature and 2 mg/L of Cd stress were measured. Results showed that there was a noticeable difference in root Cd content at the booting stage when these two rice cultivars were cultured under the moderate-temperature condition (Figure 1). Hence, YZX and XWX rice were cultivated at moderate temperature and 2 mg/L of Cd stress for 48 h to explore the gene expression alterations in response to Cd treatment in the following experiments.Figure 1
Cd accumulation difference in the roots of YZX and XWX at 4 different growth stages under different temperatures and the same Cd stress. YZX and XWX rice at different growth stages (i.e., seedling, tillering, booting, and grain-filling stages) were cultured for 48 h under different temperatures (low temperature (15–20°C), moderate temperature (22–27°C), or high temperature (30–35°C)) and 2 mg/L of Cd stress. Next, Cd concentrations in the roots of these two rice cultivars were determined.


## 3.2. Identification of Differentially Expressed Transcripts after Cd Treatment in YZX and XWX

In this study, RNA-seq technology was used to identify Cd response-related transcripts in YZX and XWX and investigate genes that might play vital roles in governing the difference in root Cd accumulation between YZX and XWX. The statistics information of clean data after filtration of RNA-seq raw data is shown in Table2. The GC content of clean data was approximately 50% (Table 2). Moreover, more than 94% of bases have a recognition accuracy of over 99.9% in the clean data (Table 2). Differential expression analysis showed that 341 transcripts were differentially expressed (206 downregulated and 135 upregulated) in the AMD vs. AMK group (Supplementary Figure 1(a) and Supplementary Table 1). Also, 75 downregulated transcripts and 86 upregulated transcripts were identified in the BMD vs. BMK group (Supplementary Figure 1(b) and Supplementary Table 2). Moreover, 18 dysregulated transcripts (i.e., BGIOSGA022091, BGIOSGA007908, BGIOSGA007302, BGIOSGA029201, BGIOSGA022060, BGIOSGA023662, BGIOSGA026675, BGIOSGA001422, BGIOSGA025476, BGIOSGA030712, BGIOSGA004724, BGIOSGA018872, BGIOSGA010590, BGIOSGA001025, BGIOSGA016478, BGIOSGA020263, BGIOSGA001117, and BGIOSGA000917) in the AMD vs. AMK group and 6 dysregulated transcripts (i.e., BGIOSGA004809, BGIOSGA030712, BGIOSGA034657, BGIOSGA012380, BGIOSGA018872, and BGIOSGA020584) in the BMD vs. BMK group were identified to be involved in inorganic ion transport and metabolism based on COG_class and KOG_class annotation analyses in Supplementary Tables 1 and 2. These transcripts might play vital roles in governing the difference in root Cd accumulation in YZX and XWX.Table 2
The statistics information of clean data.
Sample no.Clean readsClean basesGC content (%)% ≥Q30AMD-128,320,9218,469,123,35452.5194.89AMD-224,166,3847,224,223,02052.6295.04AMD-323,919,6997,147,482,54452.1695.05AMK-125,675,7837,673,764,99452.9994.99AMK-224,904,1077,449,500,70651.7595.23AMK-323,973,6977,172,841,03252.7994.77BMD-130,635,0409,159,932,34053.2896.41BMD-221,090,0586,311,270,66851.6494.84BMD-323,708,8487,088,906,52251.6794.81BMK-125,223,0397,541,950,04452.2894.88BMK-237,274,48111,138,840,83253.1396.36BMK-327,498,8598,216,506,79852.4996.44Notes. Clean reads: the total number of paired-end reads in clean data; clean bases: the total number of bases in clean data; GC content: the percentages of G and C bases in clean data; ≥Q30%: the percentage of bases with the quality value ≥ 30 (base identification accuracy ≥ 99.9%).


## 3.3. Construction of PPI Networks of Dysregulated Genes in the AMD vs. AMK and BMD vs. BMK Groups

Next, the STRING database was used to investigate the relationships of proteins encoded by the abovementioned differentially expressed genes. The interaction networks of the dysregulated genes in the AMD vs. AMK group are presented in Supplementary Figure2(a) and Supplementary Table 3. Similarly, the interaction relationships of the differentially expressed genes in the BMD vs. BMK group are shown in Supplementary Figure 2(b) and Supplementary Table 4. Additionally, the node degree (number of interacted proteins) of each protein in the PPI network of Supplementary Figures 2(a) and 2(b) was analyzed and is displayed in Supplementary Tables 5 and 6, respectively. Genes with greater node degrees might play central roles in the response to Cd. Among genes with higher node degrees, 7 genes (Os06g0196300 (OsJ_019618), Os07g0570700 (OsJ_24808), ADI1, GDCSH, HSFB2C, PEX11-4, and CLPB1) were screened out for further RT-qPCR validation. Moreover, interaction analysis of dysregulated genes in the AMD vs. AMK group revealed that there was a complex interaction among OsJ_019618, OsJ_24808, ADI1, GDCSH, and PEX11-4. For instance, OsJ_019618 could interact with OsJ_24808, ADI1, GDCSH, or PEX11-4 (Supplementary Tables 3 and 5). ADI1 could interact with OsJ_019618, GDCSH, or PEX11-4 (Supplementary Tables 3 and 5). Moreover, there was a potential interaction between OsJ_24808 and OsJ_019618 or GDCSH (Supplementary Tables 3 and 5).


## 3.4. Comparative Analysis of Cd-Responsive Genes in XWX and YZX

GO annotation analysis of the 341 differentially expressed transcripts in the AMD vs. AMK group showed that 37 transcripts were involved in the response to Cd (Supplementary Table7). Moreover, 23 of the 161 differentially expressed transcripts in the BMD vs. BMK group were identified to be implicated in Cd response based on GO annotation analysis of dysregulated transcripts in the BMD vs. BMK group (Supplementary Table 8). Among these Cd-responsive transcripts, 19 common elements were found in AMD vs. AMK and BMD vs. BMK groups (Figure 2(a) and Supplementary Table 9). These common transcripts had the same alteration trends in the AMD vs. AMK and BMD vs. BMK groups (Supplementary Table 9). And 18 or 4 transcripts were found to be differentially expressed only in the AMD vs. AMK group or BMD vs. BMK group, respectively (Figure 2(a) and Supplementary Table 9). These 22 transcripts might have crucial roles in regulating the difference in root Cd accumulation between YZX and XWX at the booting stage. Additionally, we noticed that 43 transcripts were differentially expressed only in the BMD vs. BMK group, but not notably altered in the AMD vs. AMK group (Supplementary Table 10). Also, 214 transcripts were found to be markedly upregulated or downregulated in the AMD vs. AMK group, whereas the difference in the expression of these 214 transcripts was not significant in the BMD vs. BMK group (Supplementary Table 11). The 257 transcripts in Supplementary Tables 10 and 11 are integrated into Supplementary Table 12. We believed that these transcripts might also be related to the difference in root Cd accumulation between XWX and YZX. Also, transcripts with notably different upregulated or downregulated ratios in the BMD vs. BMK and AMD vs. AMK groups might be implicated in root Cd accumulation difference between two rice cultivars, which were not analyzed in our project. KEGG enrichment analysis showed that the transcripts in Supplementary Table 12 mainly participated in the regulation of metabolic pathways, biosynthesis of secondary metabolites, carbon fixation in photosynthetic organisms, carbon metabolism, glutathione metabolism, fructose and mannose metabolism, glycolysis/gluconeogenesis, starch and sucrose metabolism, and glycerophospholipid metabolism (Supplementary Table 13 and Figure 2(b)). For instance, 23, 16, 5, 7, or 5 differentially expressed genes were significantly enriched in KEGG pathways related to metabolism, biosynthesis of secondary metabolites, carbon fixation in photosynthetic organisms, carbon metabolism, or glutathione metabolism, respectively (Supplementary Table 13). The top 20 KEGG pathways are shown in Figure 2(b).Figure 2
Comparative analysis of Cd-responsive genes in XWX and YZX. (a) Venn analysis for transcripts in Supplementary Tables7 and 8. AMD vs. AMK list: differentially expressed transcripts related to Cd response in the AMD vs. AMK group based on GO annotation analysis of dysregulated transcripts in the AMD vs. AMK group. BMD vs. BMK list: differentially expressed transcripts related to Cd response in the BMD vs. BMK group based on GO annotation analysis of dysregulated transcripts in the BMD vs. BMK group. Upper subfigure: Venn diagrams showing a number of differentially expressed transcripts related to Cd response in the AMD vs. AMK and BMD vs. BMK groups, as well as the number of overlapping genes. Middle subfigure: number of transcripts in the AMD vs. AMK and BMD vs. BMK lists. Bottom subfigure: number of specific (1) and shared (2) transcripts in the AMD vs. AMK and BMD vs. BMK lists. (b) The top 20 KEGG pathways enriched by the transcripts in Supplementary Table 12. The circle size denotes the number of transcripts in corresponding KEGG pathway terms. The circle color represents the “−log10 (P value)”.
(a)(b)


## 3.5. Expression Analysis of 10 Interested Genes by RT-qPCR Assay

Next, 10 genes including the 7 abovementioned genes with greater node degrees in the PPI networks, 1 interested gene (LAC6), and 2 genes related to heavy metal stress tolerance (WNK3 and CFBP1) (detailed information is shown in the Discussion section) were selected for further RT-qPCR validation. RT-qPCR results showed that the expression levels of ADI1, CFBP1, GDCSH, LAC6, Os06g0196300, Os07g0570700, PEX11-4, and WNK3 were notably downregulated in both AMD vs. AMK and BMD vs. BMK groups (Figure3(a), 3(b), 3(d), and 3(f)–3(j)). CLPB1 and HSFB2C were highly expressed in the AMD vs. AMK group, but not notably changed in the BMD vs. BMK group (Figures 3(c) and 3(e)). Comparative analysis of RT-qPCR and RNA-seq outcomes disclosed that the alteration trends of these 10 genes were consistent in response to Cd treatment in YZX and XWX rice cultivars, although some differences were not statistically significant (Figures 3(a)–3(j)).Figure 3
Comparative analysis of RT-qPCR and RNA-seq outcomes of expression patterns of 10 interested genes (i.e., PEX11-4 (a), WNK3 (b), CLPB1 (c), GDCSH (d), HSFB2C (e), LAC6 (f), Os06g0196300 (g), Os07g0570700 (h), ADI1 (i), and CFBP1 (j)) in roots of XWX and YZX.
(a)(b)(c)(d)(e)(f)(g)(h)(i)(j)


## 3.6. Interaction and Enrichment Analysis of These 10 Genes

The relationships of these 10 genes were further deciphered by PPI analysis via the STRING database. Results suggested that 7 gene-encoded proteins (ADI1, CFBP1, PEX11-4, OsJ_019618, OsJ_24808, GDCSH, and CLPB1) could form a potential regulatory network (Figure4). The interaction analysis of these 10 proteins is shown in sheet 1 of Supplementary Table 14. Moreover, GO enrichment analysis revealed that 7 genes (LAC6, CFBP1, OsJ_019618, WNK3, OsJ_24808, ADI1, and GDCSH) were enriched in the metabolic process (Supplementary Table 14, sheet 2). Additionally, 6 genes (LAC6, CFBP1, CLPB1, WNK3, OsJ_24808, and ADI1) had a potential ion binding activity (Supplementary Table 14, sheet 3). Furthermore, annotation analysis revealed that 3 genes (LAC6, OsJ_019618, and PEX11-4) were involved in the regulation of the response to oxidative stress (Supplementary Table 14, sheet 4).Figure 4
Interaction network of 10 interested genes. The yellow line represents that the potential interaction between proteins was obtained by text mining. The black line denotes that the potential interaction between proteins was identified by coexpression analysis. The green line means that the potential interaction between proteins was predicted based on gene neighborhood relationships. PPI analysis suggested that PEX11-4 could interact with CFBP1 (FBP1), ADI1, and Os06g0196300 (OsJ_019618). There was a potential interaction between CFBP1 (FBP1) and PEX11-4/ADI1/Os06g0196300 (OsJ_019618)/GDCSH. GDCSH coud interact with CFBP1 (FBP1)/ADI1/Os06g0196300 (OsJ_019618)/Os07g0570700 (OsJ_24808). Os07g0570700 (OsJ_24808) could interact with GDCSH/Os06g0196300 (OsJ_019618). Os06g0196300 (OsJ_019618) could interact with ADI1/PEX11-4/CFBP1 (FBP1)/GDCSH/Os07g0570700 (OsJ_24808)/CLPB1. ADI1 could interact with PEX11-4/CFBP1 (FBP1)/GDCSH/Os06g0196300 (OsJ_019618). CLPB1 could interact with Os06g0196300 (OsJ_019618).


## 3.7. Correlation Analysis of These 10 Genes and Cd Response

As shown in Supplementary Table5, OsJ_019618, OsJ_24808, ADI1, GDCSH, PEX11-4, CLPB1, and HSFB2C could interact with 59, 44, 36, 29, 17, 6, and 6 genes in the interaction network of dysregulated genes in the AMD vs. AMK group, respectively. And CLPB1, PEX11-4, and HSFB2C could interact with 6, 5, and 3 genes in the BMD vs. BMK group, respectively (Supplementary Table 6). To explore the association of these genes and Cd response, Venn analysis of genes that could interact with the 10 abovementioned genes and GO-annotated Cd-related genes in Supplementary Tables 7 and 8 was performed. Results showed that OsJ_019618 could interact with 4 genes related to the Cd response (i.e., APX2, HCF136, GLN2, and ALDP). ADI1 could interact with 2 Cd-related genes (i.e., ALDP and GLN2). OsJ_24808, GDCSH, or CLPB1 could interact with Cd-responsive gene TPI, GLN2, or HSP17.4, respectively. These data further suggested the close association of these genes and Cd response.


## 4. Discussion

Rice is a big reservoir of potentially toxic heavy metals in paddy soil-rice systems [22, 23]. The contamination of heavy metals in soil not only influences rice growth and grain yields but also seriously threatens the health of animals and humans that consume rice [22–24]. In this study, genes related to Cd response and root Cd content difference between YZX and XWX were further examined in the roots of these two rice cultivars with different grain Cd accumulation.Both RT-qPCR and RNA-seq outcomes suggested that 10 genes (i.e., PEX11-4, CLPB1, HSFB2C, LAC6, CFBP1, OsJ_019618, WNK3, OsJ_24808, ADI1, and GDCSH) might play vital roles in the response to Cd in both YZX and XWX. Moreover, interaction analysis suggested that 7 genes (ADI1, CFBP1, PEX11-4, OsJ_019618, OsJ_24808, GDCSH, and CLPB1) could form a potential regulatory network in rice. Additionally, OsJ_019618 could interact with 4 GO-annotated genes related to Cd response (i.e., APX2, HCF136, ALDP, and GLN2). Moreover, there was a potential interaction between ADI1 and GO-annotated Cd-related gene ALDP or GLN2. OsJ_24808, CLPB1, or GDCSH could interact with GO-annotated Cd-responsive gene TPI, HSP17.4, or GLN2, respectively. Furthermore, APX2 [25, 26], HCF136 [27], ALDP [28], TPI [29, 30], and HSP17.4 [31] have been found to be implicated in Cd response. ADI1 has also been reported to be related to Cd stress in Agaricus brasiliensis [32]. Cd exposure could lead to the increase of CLPB expression level and delay of cell division in Escherichia coli [33, 34]. Additionally, our analysis revealed that 6 genes (LAC6, CFBP1, CLPB1, WNK3, OsJ_24808, and ADI1) had a potential ion binding activity. Probable serine/threonine-protein kinase WNK3 belongs to the WNK family, which plays vital roles in plant growth and development and stress responses [35, 36]. A recent study showed that WNK9 (a member of the WNK family) overexpression could improve arsenite (a heavy mental) tolerance in transgenic Arabidopsis [37]. Furthermore, annotation analysis revealed that 3 genes (LAC6, OsJ_019618, and PEX11-4) were related to oxidative stress. Oxidative stress has been found to be closely linked with Cd tolerance and Cd-induced injury and toxicity [38–40]. For instance, molybdenum could relieve Cd-induced toxicity and potentiate Cd tolerance by restricting Cd uptake, reducing oxidative stress, and improving antioxidant defense responses [41]. These data suggested that these genes might function as crucial players in the response to Cd stress.Cao et al. showed that 706 transcripts were specially expressed in the roots of the mutantIndica rice (Icd1 group), which had an extremely low Cd accumulation in root and shoot under Cd exposure [42]. Moreover, 987 transcripts were specially expressed in the roots of wild-type rice (WT group) under Cd exposure [42]. These transcripts might be implicated in the difference in Cd accumulation in wild-type and mutant-type Icd1 rice. Combined with our data and Cao’s outcomes, we found that 4 transcripts (BGIOSGA004865 (CFBP1), BGIOSGA017856 (Os05g0480200), BGIOSGA031509 (GSTU6), and BGIOSGA033278 (GDCSH)) that were specially expressed in the Icd1 group were differentially expressed in the AMD vs. AMK group. Also, 11 transcripts (BGIOSGA000239 (SPS1), BGIOSGA000945 (LEA14), BGIOSGA001233 (ABCG36), BGIOSGA005931 (YSL2), BGIOSGA006919 (SAT1), BGIOSGA007802 (CKX6), BGIOSGA011045 (HSP17.4), BGIOSGA014421 (AOX1B), BGIOSGA021670 (OPR1), BGIOSGA033315 (GSTU6), and BGIOSGA033325 (GSTU6)) or 6 transcripts (BGIOSGA000239 (SPS1), BGIOSGA006919 (SAT1), BGIOSGA009154 (HSP18.6), BGIOSGA011045 (HSP17.4), BGIOSGA021670 (OPR1), and BGIOSGA033325 (GSTU6)) that were specially expressed in the WT group were differentially expressed in the AMD vs. AMK or BMD vs. BMK group, respectively. Among these genes, 2 genes (GSTU6 and HSP17.4) were annotated to be involved in Cd response by the GO database. Moreover, GSTU6 has been found to be notably upregulated in roots of rice under Cd stress [43]. HSP17.4 expression was markedly increased in Arabidopsis thaliana roots after Cd exposure [31]. SPS1 was notably downregulated in response to Cd exposure in Chlorella sorokiniana [44]. ABCG36 (PDR9) expression was markedly increased in the roots of rice after a short exposure to Cd [45, 46]. ABCG36 knockout promoted root Cd accumulation and potentiated Cd sensitivity in rice [45]. The YSL2 expression level was noticeably increased in the root and shoot of rice in response to Cd stress and YSL2 might be related to Cd uptake transport [47, 48]. These transcripts might play vital roles in the response to Cd stress in rice.


## 5. Conclusion

Taken together, a host of Cd response-related transcripts were identified in the roots of YZX and XWX at the booting stage. Some transcripts that might be associated with the difference in the root Cd-responsive ability of these two rice cultivars were screened out. Moreover, our data suggested that ADI1, CFBP1, PEX11-4, OsJ_019618, OsJ_24808, GDCSH, CLPB1, and WNK3 might play vital roles in regulating Cd response and Cd tolerance. An in-depth insight into the genetic and molecular mechanisms underlying root Cd stress response difference in YZX and XWX at the booting stage might contribute to the development of new strategies that can improve rice yield and reduce Cd harm to rice and humans.



---

*Source: 2166775-2021-12-13.xml*