Dissecting Protein-Protein Interaction Networks of Arabidopsis thaliana and Arabidopsis halleri to Get Insights into Heavy Metal Tolerance Strategies †

: To gain insight into two different plant strategies (hyperaccumulation and phytostabilization) for managing heavy metals, we conducted a network based functional enrichment analysis. Protein-protein interactions of A. halleri root and shoot were derived by weighted gene co-expression analysis. While, in the case of A. thaliana , protein-protein interactions characterizing the organs were derived from STRING database based on genes known to be expressed in root and shoot. Protein-protein interaction clusters of root and shoot networks of both species were analyzed to identify enriched pathways. Thus, we provide a ﬁrst clear analysis of the biological peculiarities of different organs of both species.


Introduction
Phytoremediation is one of the auspicious strategies to manage heavy-metal polluted soils [1]. Several plants act differently to accommodate this issue rather by uptake (phytoextraction), emit in the atmosphere (phytovolatilization), or stabilize heavy metals in the root system (phytostabilization) [2]. Among these different processes, phytostabilization is of particular interest to be a low-cost and effective strategy. Moreover, this process is characteristic of a wide range of plants which include also commercial species and model organisms such as Arabidopsis thaliana [3]. On the other hand, some plants developed a rare adaptation to extract heavy metals from the soil and hyperaccumulate them into the shoot, such in the case as Arabidopsis halleri [4]. Within the genus Arabidopsis, A. thaliana and A. halleri possesses an unique global distributions, mating systems, life histories and adaptations providing opportunities to study them under different growth conditions. The recently long-scaffold assembly of the A. halleri permits now to identify long-range patterns of polymorphism and diversity and perform further genotype-phenotype association studies in different fields [5]. The major area of research in A. halleri focused on the study of heavy metal tolerance and hyperaccumulation, a constitutive phenotype in all tested genotypes. However, up to now, molecular mechanisms underpinning this ability are not completely understood; thus, in the present study, the prediction of all protein-protein interactions (PPIs), followed by network based functional enrichment analysis, were used to identify all pathways characterizing A. halleri, followed by comparing the uniqueness and similarities with A. thaliana.

Methodology
A weighted gene co-expression network analysis (WGCNA) was applied on RNA-seq data of A. halleri root and shoot [6] to identify and describe the PPI for the first time by WGCNA, a free accessible R package [7]. A soft threshold power was identified for root and shoot, separately, and used it to construct adjacency matrix based on scale-free topology. The adjacency matrix was converted to topological overlap matrix (TOM) and related dissimilarity of TOM (1-TOM) was computed to classify genes, having same expression pattern, into same module. Clustering height cut-off was set to 0.25 in order to merge likely modules with minimum module size of 30.
For A. thaliana, taking advantage of already available data, a list of genes expressed in the root and shoot was identified by reviewing the recent published papers, and organspecific PPI was reconstructed and, successively, clustered according to their topological position into the networks. The gene lists were used to filter a high-confidence set of interactions of A. thaliana derived from STRING database [8]. The resulting network was subjected to MCODE plugin of Cytoscape to identify the highly connected regions in the PPI network which represent molecular clusters (modules) [9].
After identifying and filtering the clusters/modules in the two species, we performed KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis by using g:Profiler to have a view of biological pathways in common and dissimilar and involvement of genes [10].

Network Construction of A. halleri
The WGCNA was carried out on RNA-seq data of A. halleri root and shoot. A soft threshold power of β = 14 for root, and β = 8 in case of shoot was used to construct adjacency matrix based on scale-free topology (R2 = 0.97) and, successively, all genes having similar expression patterns were categorized into 14 modules in root and shoot ( Figure 1). pathways characterizing A. halleri, followed by comparing the uniqueness and similarities with A. thaliana.

Methodology
A weighted gene co-expression network analysis (WGCNA) was applied on RNAseq data of A. halleri root and shoot [6] to identify and describe the PPI for the first time by WGCNA, a free accessible R package [7]. A soft threshold power was identified for root and shoot, separately, and used it to construct adjacency matrix based on scale-free topology. The adjacency matrix was converted to topological overlap matrix (TOM) and related dissimilarity of TOM (1-TOM) was computed to classify genes, having same expression pattern, into same module. Clustering height cut-off was set to 0.25 in order to merge likely modules with minimum module size of 30.
For A. thaliana, taking advantage of already available data, a list of genes expressed in the root and shoot was identified by reviewing the recent published papers, and organspecific PPI was reconstructed and, successively, clustered according to their topological position into the networks. The gene lists were used to filter a high-confidence set of interactions of A. thaliana derived from STRING database [8]. The resulting network was subjected to MCODE plugin of Cytoscape to identify the highly connected regions in the PPI network which represent molecular clusters (modules) [9].
After identifying and filtering the clusters/modules in the two species, we performed KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis by using g:Profiler to have a view of biological pathways in common and dissimilar and involvement of genes [10].

Network Construction of A. halleri
The WGCNA was carried out on RNA-seq data of A. halleri root and shoot. A soft threshold power of β = 14 for root, and β = 8 in case of shoot was used to construct adjacency matrix based on scale-free topology (R2 = 0.97) and, successively, all genes having similar expression patterns were categorized into 14 modules in root and shoot ( Figure 1).

Network Reconstruction of A. thaliana
Among all STRING interactions available for A. thaliana, we considered only those with a score higher than 700 composing a network made of a high-confidence set of interactions and only if both of its node IDs were present in the list of genes expressed in the two organs [11]. The modules (clusters) within each organ-specific networks were

Network Reconstruction of A. thaliana
Among all STRING interactions available for A. thaliana, we considered only those with a score higher than 700 composing a network made of a high-confidence set of interactions and only if both of its node IDs were present in the list of genes expressed in the two organs [11]. The modules (clusters) within each organ-specific networks were identified by applying a molecular complex detection algorithm available as MCODE plugin of Cytoscape, with the following criteria: degree cut-off equal to 2, node score cut-off equal to 0.2, K-core score equal to 2 and maximum depth from the seed equal to 100. A final score was calculated for the i-th module according to the formula: where D i is the density of the module and N i the number of nodes of the module. Modules were ranked and retained for following analysis if their final score was higher than 10. By this procedure we identified 15 modules including 1023 genes from root network and 18 modules including 1804 genes form shoot.

KEGG Enrichment Analysis
With the purpose of getting insights of the biological pathways characterizing the biology of root and shoot of A. halleri and A. thaliana, KEGG pathway enrichment analysis was performed and identified the top 3 KEGG pathways significantly enriched (p < 0.05) in each module of root and shoot in A. halleri and A. thaliana (Figure 2). identified by applying a molecular complex detection algorithm available as MCODE plugin of Cytoscape, with the following criteria: degree cut-off equal to 2, node score cutoff equal to 0.2, K-core score equal to 2 and maximum depth from the seed equal to 100. A final score was calculated for the i-th module according to the formula: where D i is the density of the module and N i the number of nodes of the module. Modules were ranked and retained for following analysis if their final score was higher than 10. By this procedure we identified 15 modules including 1023 genes from root network and 18 modules including 1804 genes form shoot.

KEGG Enrichment Analysis
With the purpose of getting insights of the biological pathways characterizing the biology of root and shoot of A. halleri and A. thaliana, KEGG pathway enrichment analysis was performed and identified the top 3 KEGG pathways significantly enriched (p < 0.05) in each module of root and shoot in A. halleri and A. thaliana (Figure 2). In case of A. halleri, the most enriched pathway in both root and shoot was "photosynthesis". While, in case of A. thaliana, the pathway related to "ribosome" was In case of A. halleri, the most enriched pathway in both root and shoot was "photosynthesis". While, in case of A. thaliana, the pathway related to "ribosome" was most enriched in both organs. Figure 3 is showing a Venn diagram representing the number of common and unique pathways identified in the root and shoot of both plants. In case of root, there are 8 unique pathways for A. halleri (24.2% of the total pathways identified), whereas for A. thaliana, there were 16 unique pathways (48.5% of total). On the other hand, in shoot, there were 5 unique pathways in A. halleri (14.7% of the total) and 24 pathways in case of A. thaliana (70.6% of the total).
Biol. Life Sci. Forum 2021, 11, 13 4 of 6 most enriched in both organs. Figure 3 is showing a Venn diagram representing the number of common and unique pathways identified in the root and shoot of both plants.  Glucosinolate biosynthesis was one of the enriched unique pathways identified in the in root of A. halleri. Glucosinolates comprise of diverse range of organic compounds involved in interactions between plant and insects [12]. Few studies reported the role of glucosinolates in anti-herbivory in A. thaliana and A. halleri [13,14]. Cadmium exposure to A. halleri results in the increased accumulation of glucosinolates in leaves and the result was opposite in case of A. thaliana with the decreased deposition of some glucosinolates in leaves and roots with Cd exposure [13,15]. These results can be relatable to our results where glucosinolate pathway was found to be more enriched in root of A. halleri. On the other hand, in A. thaliana root, citrate cycle (TCA cycle) was found to be among the unique pathways. TCA cycle is necessary to maintain the normal growth and development of plant under stress environment [16]. A study revealed that Cd exposure caused an increased accumulation of enzymes involved in TCA cycle, including citrate synthase (At2g44350) [17], which is also identified to be involved in this pathway in our study. The significance of this pathway in A. thaliana is well defined also in literature, as knockout of one of its components resulted in dwarfing phenotype accompanying with more prone to oxidative stress [18]. Till date, there is no prominent study regarding the component of this pathway in A. halleri. However, a study found more citrate (produced in TCA cycle) in the roots of A. halleri grown in contaminated soil with Zn [5,19]. But here we cannot conclude based on a single result.
Phenylpropanoids biosynthesis, pointed out as unique in case of A. halleri shoot, is a pathway which consists of a sequence of enzyme regulated reactions leading to different aromatic end products [20]. A study summarized the involvement of this pathway in increased tolerance against several abiotic stresses including heavy metals [21]. However, there is a lack of knowledge regarding the relation between metal stress and phenylpropanoid biosynthesis pathway in A. halleri and A. thaliana. But it should be considered that there is a chance of possible crosstalk between phenylpropanoid and glucosinolate metabolic pathways [22,23]. Based on this view, we can expect the same trend of phenylpropanoid metabolic pathway in accordance with the previous mentioned review in terms of glucosinolate metabolic pathway. For A. thaliana, pentose phosphate pathway and SNARE interactions were among the most enriched unique pathways in shoot. The pentose phosphate pathway is mandatory for metabolism as it produces components essential for nucleotide synthesis and for aromatic amino acids [24]. However, according to our knowledge, no study has been reported regarding the role of this pathway in metal stress in both plants. On the other side, several studies mentioned Glucosinolate biosynthesis was one of the enriched unique pathways identified in the in root of A. halleri. Glucosinolates comprise of diverse range of organic compounds involved in interactions between plant and insects [12]. Few studies reported the role of glucosinolates in anti-herbivory in A. thaliana and A. halleri [13,14]. Cadmium exposure to A. halleri results in the increased accumulation of glucosinolates in leaves and the result was opposite in case of A. thaliana with the decreased deposition of some glucosinolates in leaves and roots with Cd exposure [13,15]. These results can be relatable to our results where glucosinolate pathway was found to be more enriched in root of A. halleri. On the other hand, in A. thaliana root, citrate cycle (TCA cycle) was found to be among the unique pathways. TCA cycle is necessary to maintain the normal growth and development of plant under stress environment [16]. A study revealed that Cd exposure caused an increased accumulation of enzymes involved in TCA cycle, including citrate synthase (At2g44350) [17], which is also identified to be involved in this pathway in our study. The significance of this pathway in A. thaliana is well defined also in literature, as knockout of one of its components resulted in dwarfing phenotype accompanying with more prone to oxidative stress [18]. Till date, there is no prominent study regarding the component of this pathway in A. halleri. However, a study found more citrate (produced in TCA cycle) in the roots of A. halleri grown in contaminated soil with Zn [5,19]. But here we cannot conclude based on a single result.
Phenylpropanoids biosynthesis, pointed out as unique in case of A. halleri shoot, is a pathway which consists of a sequence of enzyme regulated reactions leading to different aromatic end products [20]. A study summarized the involvement of this pathway in increased tolerance against several abiotic stresses including heavy metals [21]. However, there is a lack of knowledge regarding the relation between metal stress and phenylpropanoid biosynthesis pathway in A. halleri and A. thaliana. But it should be considered that there is a chance of possible crosstalk between phenylpropanoid and glucosinolate metabolic pathways [22,23]. Based on this view, we can expect the same trend of phenylpropanoid metabolic pathway in accordance with the previous mentioned review in terms of glucosinolate metabolic pathway. For A. thaliana, pentose phosphate pathway and SNARE interactions were among the most enriched unique pathways in shoot. The pentose phosphate pathway is mandatory for metabolism as it produces components essential for nucleotide synthesis and for aromatic amino acids [24]. However, according to our knowledge, no study has been reported regarding the role of this pathway in metal stress in both plants. On the other side, several studies mentioned the role of SNAREs in abiotic tolerance in plants. The decrease in expression of v-SNAREs proteins in A. thaliana resulted in improved salt tolerance, leading to efficient functioning of tonoplast and vacuole [25].
This aspect can also be considered to regulate redox activity in relation to heavy metal stress, as the heavy metals have the ability to generate reactive oxygen species [26]. However, there is a lack in knowledge of these proteins in A. halleri.

Conclusions
Our study presented, for the first time, the prediction of PPI in A. halleri and, followed by the functional enrichment analysis, we provide new insights into the biology of this recently sequenced organism. In addition, a simple methodology was established to gain knowledge about similar or contrasting biological behaviors of A. halleri and A. thaliana.
A. halleri is a close relative of A. thaliana, but the main mechanisms for tolerating metal stress may be specific to one or the other species and therefore may be based on the possession of some different and particular biological processes. Here, we pointed out few organ-specific biological pathways identified as unique to A. halleri and/or A. thaliana, which can be suggested as possible targets to do future research. Moreover, this study will also help the researchers in order to pinpoint the gaps in research and to explore more of these two Arabidopsis species.

Data Availability Statement:
The data underlying this article will be shared on reasonable request 236 from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.