BACKGROUND
as an essential constituent of many important compounds that are required for the development and growth of plants, phosphorus  actively participates in photosynthesis, respiration, carbohydrate metabolism, energy transduction and other processes  <cit> .
in many agricultural systems, however, the concentration of absorbable phosphorus in soil is insufficient to ensure plant productivity  <cit> .
as one of the most important grain and industrial crops, soybean  is a vital source of protein for human and animal food.
however, soybean growth and yield are limited by varieties of edaphic factors, especially low phosphorus availability in soils  <cit> .
under low-p stress, soybean plants become dwarfed, leaf area is reduced, and necrotic spots might appear on lower leaves  <cit> .
also, flowering time is increased and the number of pods could decrease during the fruiting period  <cit> .
due to these p-mediated phenotypic effects on soybeans, phosphate fertilizers are used to overcome undesirable potential yield outcomes.
large quantities of phosphate fertilizer may cause environmental pollution as it is easily immobilized in the soil to form chelates or precipitates  <cit> .
in order to alleviate the effects of phosphorus deficiency on soybean growth, it is particularly important for us to search for p-efficient soybean materials and investigate high-phosphorus-efficient genes.
recently, four major groups of genes induced by low-p stress were reported  <cit> .
one group was involved in sensing low-p signals, including local phosphate sensing that was proposed to occur in the primary root in arabidopsis thaliana  <cit> , and long-distance phosphate signaling which depends on the activity of transcription factors, such as phr <dig>  phl <dig>  <cit> .
members of another group mediated phosphate distribution by involving an internal regulation of phosphorus homeostasis, including spx family and pho <dig> that contains spx domain in arabidopsis  <cit> .
the remaining two groups of genes mainly focus on the acquisition and transportation of phosphate.
for example, gmacp <dig> encodes an acid phosphatase, the overexpression of which in roots could increase the absorption of phosphate in the nutrient solution  <cit> .
the pht <dig> family consisting of  <dig> pht <dig> genes in soybean encoded phosphate transporters was demonstrated to be involved in not only directing phosphate uptake in roots, but also transporting phosphorus as phosphate  <cit> .
generally, the acquisition and transportation of phosphate greatly depends on the root system as root morphology changes under p-deficient conditions to sustain plant growth.
this includes an increase in number and length of lateral roots in maize   <cit> , longer and greater number of root hairs in low-p tolerant accessions of barley  and arabidopsis  <cit> , and relatively longer primary roots in rice  when a low-phosphorus-related transcription factor osmyb2p- <dig> is over expressed  <cit> .
meanwhile, qtl  mapping associated with soybean phosphorus efficiency have made processes.
qtls related to seed phosphorus content  <cit> , phosphorus efficiency relevant root traits  <cit> , and soybean tolerance to low phosphorus stress based on flower and pod abscission rate  <cit>  have been identified.
furthermore, screening of soybean p-efficient materials identified p-efficient genotypes, such as hp <dig>  hp <dig>  huaxai  <dig>  huaxia  <dig> etc.
<cit> .
as an abiotic stress, low-p conditions may induce a series of complex responses.
it was demonstrated that phytohormone signaling pathways responded strongly to phosphate deprivation  <cit> .
the expression of a large set of genes associated with plant hormone signaling either induced or repressed root growth in arabidopsis and zea mays when suffering low-p conditions  <cit> .
in addition, changes of the secondary metabolites to low-p stress were shown with an accumulation of anthocyanin under long-term phosphorus deficiency  <cit> .
similar results were observed in zea mays with genes related to the phenylpropanoid pathway identified as up-regulated or down-regulated under p starvation  <cit> .
a large set of genes and signaling pathways that were associated with low-p stress had been uncovered using the technique of gene expression profiles and transcriptome analyses in arabidopsis, rice, and zea mays .
recently, zeng et al.
identified phosphate-deficiency-responsive genes in soybean  roots by high-throughput sequencing  <cit> .
selection of representative materials and consideration of different tissues will further reveal genes and pathways associated with low-p stress.
here, we selected a low-p-tolerant accession chundou  and a low-p-sensitive accession yunhefengwodou  from  <dig> soybean accessions through hydroponic experiments.
then microarray chips were performed on the roots and leaves of these two accessions.
by analyzing the gene expression data, we identified  <dig> candidate genes and three common pathways that were induced by low-p stress.
our study provides candidate loci for functional identification of high-phosphorus-efficient genes which may be of great significance for cultivating p-efficient soybean accessions.
RESULTS
identification of low-p-tolerant accession cd and low-p-sensitive accession yh
with hydroponic experiments lasting for 10 days under different p treated conditions of  <dig> soybean accessions  <cit> , we found that the soybean accession cd and yh showed similar root morphology under + p treated condition , while showing obvious morphology differences under -p treated condition .
we then identified accession cd as a low-p-tolerant accession, and accession yh as a low-p-sensitive accession.
then, we carried out a gradient experiment using full, half and quarter hoagland’s solution and found that the + p and -p treated soybean plants showed most obvious morphology difference when grown in half hoagland.
hence half hoagland with  <dig>  mmol/l p and  <dig>  mmol/l p were considered as normal and low p treated conditions, respectively.
the phenotypic traits of cd and yh are shown in additional file 1: figure s <dig>  primary root length was longer in cd than that in yh.
the number of root hairs and lateral roots of cd were much more than those in yh.
microarray experiments were performed on the roots and leaves of cd and yh.
significantly differentially expressed genes in different p treated conditions
under low p condition,  <dig> and  <dig> up-regulated differently expressed genes  were found in the roots and leaves of low-p-tolerant accession cd respectively; while no up-regulated degs were found in the roots and leaves of low-p-sensitive accession yh compared with low-p treatment.
for down-regulated degs, only  <dig> degs could be found in the roots of cd,  <dig> degs in the roots, and  <dig> degs in the leaves of yh .table  <dig> the number of degs in different p treated conditions
generally, the number of degs in cd was more than that in yh, not only in the roots but also in the leaves , which suggested that cd could regulate a series of genes to cope with the low-p stress.
we found only one common deg, glyma20g <dig> , in both the roots of cd and yh, with the fold change  <dig>  and  <dig>  down-regulated in low-p-tolerant accession cd and low-p-sensitive accession yh, respectively.
the number of degs in the roots was more than that in the leaves.
interestingly, one common deg, glyma01g <dig>  named as gmmmp <dig>  which encodes a matrix metalloproteinase, was found both in the roots and leaves of cd.
the fold change is  <dig>  and  <dig> , respectively.
previous studies demonstrated that it was associated with pathogenic infections in soybean  <cit> .
a total of  <dig> non-redundant degs were found between treated and control conditions in the roots and leaves of both soybean accessions used in this experiment.
the hierarchical clustering analysis classified the tolerant and sensitive soybean accessions into two clusters; samples cluster and genes cluster.
the compactness of the clusters of the three replicates confirmed the high reliability of our experimental design .fig.
<dig> hierarchical clustering of  <dig> differentially expressed genes in different p-treated conditions.
a <dig>  a <dig> and a <dig>  three biological replicates of control groups ; a <dig>  a <dig> and a <dig>  three biological replicates of treatment groups .
red and green indicate high and low expression levels, respectively
to evaluate the potential functions of the degs between low-p stress and normal conditions in both soybean accessions, gene ontology  categories were applied for low-phosphorus-induced genes.
the three go terms were categorized as biological processes, molecular functions, and cellular components.
we found that degs were significantly associated with response to wounding, response to jasmonic acid stimulus, oxidoreductase activity, and cell wall, whether they were up-regulated genes or down-regulated genes.
furthermore, we identified pathways that were associated with low-p induced genes with the use of kegg .
overall,  <dig> pathways with p values less than  <dig>  were motivated under low-p stress , and most pathways were related to metabolic process.
notably, the most significant  <dig> pathways that were enriched included phenylpropanoid biosynthesis, methane metabolism and phenylalanine metabolism, all related to the metabolism of secondary metabolites.
significantly differentially expressed genes between different soybean accessions
accessions cd and yh showed strong differences in the ability to suffer low-p stress .
there were much more degs in cd  than in yh  .
to further reveal the degs related with low-p condition between cd and yh, we compared the degs between soybean accessions cd and yh .table  <dig> the number of degs in cd/yh and roots/leaves
low p ↑a
normal p → c
normal p ↓ d
a significantly up-regulated expression pattern under low-p stress;
b significantly up-regulated expression pattern under normal p condition;
c the expression patterns were not significantly changed under normal p condition;
d significantly down-regulated expression pattern under normal p condition;
e the number in parentheses indicates the number of genes identified by our criteria
degs induced by low-p conditions were classified as three categories according to their expression patterns in normal conditions.
for  <dig> up-regulated genes in roots between cd and yh in low-p stress , the three categories were as follows: 1)  <dig> genes were up-regulated in cd/yh in normal conditions, which is consistent with the expression pattern in low-p stress; 2)  <dig> genes were down-regulated in normal conditions, which showed the opposite expression patterns between normal and low-p conditions; and 3)  <dig> genes showed no expression difference in normal condition, with the expression difference less than two fold between cd and yh.
by the criterion that genes with two fold expression change between experimental group and control group were identified as degs, the expression level in cd was more than two fold higher than that in yh for the  <dig> up-regulated genes .
among the  <dig> up-regulated genes, there were  <dig> genes with no expression difference in normal conditions, the expression difference of which between cd and yh was less than two fold.
to make the expression difference of cd/yh between low-p stress and normal conditions higher than two fold, we chose those genes with larger than 4-fold expression difference between cd and yh under low-p conditions as candidates.
by this criterion,  <dig> genes with larger than 4-fold expression difference between cd and yh under low-p condition from the  <dig> genes were chosen as candidates.
hence, we got  <dig>   <dig>  and  <dig> degs from  <dig>   <dig>  and  <dig> degs, respectively .
combined with the genes with opposite expression patterns, they were considered as marked degs between cd and yh.
therefore,  <dig> marked degs were obtained.
among these genes, one gene was up-regulated more than  <dig> fold between cd/yh under low-p conditions in both roots and leaves.
thus, there were  <dig> non-redundant marked degs detected in total.
among these marked degs, the majority  were up-regulated in cd/yh under low-p stress.
additionally, there were  <dig> and  <dig> genes expressed differently between cd and yh in roots and leaves, respectively.
this result suggested more genes were activated in roots than in leaves under low-p stress during seedling stage.
we performed the enrichment of go and kegg analysis for these  <dig> marked degs.
go analysis showed that these genes were mainly involved in secondary metabolic processes, for example flavonoid biosynthesis.
they were also involved in response to hormone stimuli , and defense and abiotic stimuli .
in addition, they were involved in electron carrier activity and peroxidase.
kegg analysis revealed the pathways enriched were mostly associated with the metabolism of secondary metabolites, including phenylpropanoid biosynthesis, methane metabolism, phenylalanine metabolism, flavonoid biosynthesis and so on .
the  <dig> degs in different p-treated conditions were considered as genes in response to low-p stress.
the  <dig> marked degs in different soybean accessions were not only induced by low-p stress but also by materials-variation dependent upon the response to low-p conditions.
thus we compared  <dig> degs and  <dig> marked genes, finding  <dig> overlapping genes that were considered as active genes induced by low-p stress.
further study revealed most of the common genes focused on roots, which was consistent with a greater number of degs in roots than leaves in different p conditions.
interestingly, the  <dig> genes that showed opposite expression patterns in roots all belonged to the overlapping genes .
despite two genes  that could not find corresponding symbol and one gene  that could not find functional annotation, seven genes  with functional annotations were shown to be associated with biotic and abiotic stresses.
likewise, go and kegg enrichment analyses were applied to evaluate the potential functions of the  <dig> overlapping degs.
degs significantly enriched in go terms of biological processes such as oxidation reduction and response to jasmonic acid stimulus.
molecular functions such as oxidoreductase activity, and cellular components such as the cell wall, were also enriched.
meanwhile, sixteen pathways were enriched based on  <dig>  significance level.
similar to the results of the  <dig> marked genes, the pathways of phenylpropanoid biosynthesis, methane metabolism, phenylalanine metabolism, and flavonoid biosynthesis were also detected .
significantly differentially expressed genes between roots and leaves
in our study, under low-p treatment, more degs were detected in roots than leaves, suggesting there were more genes in response to low-p stress during seedling stage in roots.
degs between roots and leaves were found .
degs induced by low-p conditions were classified into three categories according to their expression patterns in normal conditions.
for  <dig> up-regulated genes between roots and leaves in cd under low-p stress, three categories of degs were  <dig> up-regulated genes, one down-regulated gene, and  <dig> genes showing no different expression in roots/leaves in normal conditions.
among the  <dig> genes, there were  <dig> genes that showed larger than 4-fold expression difference between roots and leaves under low-p conditions.
additionally, we found  <dig>   <dig>  and  <dig> genes showing larger than 4-fold expression differences under low-p conditions from the  <dig>   <dig>  and  <dig> genes, respectively .
with one opposite-expressed gene, there were  <dig> marked degs between roots and leaves.
we performed the enrichment of go and kegg analysis for these  <dig> marked degs.
go analysis showed these genes enriched in a group of biological processes, for example, oxidation reduction, response to oxidative stress, and response to abiotic stimulus.
kegg analysis revealed most of these genes enriched in methane metabolism, phenylpropanoid biosynthesis, and phenylalanine metabolism .
the  <dig> marked degs between roots and leaves were related with low-p stress.
there were  <dig> degs in roots and leaves considered to be caused by low-p stress.
thus, we compared  <dig> degs and  <dig> marked genes and found  <dig> overlapping genes.
interestingly,  <dig> overlapping genes were up-regulated in cd between roots and leaves, consistent with the results that more genes were induced in cd than yh under low-p stress.
go enrichment analysis showed the  <dig> overlapping genes enriched in biological processes such as response to cadmium ion, defense, jasmonic acid stimulus .
in the molecular function term and cellular component, these genes were significantly enriched in transcription regulator activity and in the nucleus  .
furthermore, we identified pathways including methane metabolism and flavonoid biosynthesis .
common degs among treatments, materials, and tissues
the  <dig> marked degs between cd and yh under low-p stress were identified as degs between materials which were associated with low-p stress, while the  <dig> marked degs between roots and leaves in low-p conditions were regarded as degs between tissues that were related to low-p stress.
we compared the  <dig>   <dig> marked degs with the  <dig> degs between different p treatments and  <dig> common degs were found to be in all three classes.
that is to say, the  <dig> degs occurred not only between treatments, but also between materials and tissues under low-p stress, which suggested they were active genes under low-p conditions.
further kegg analyses enriched in three general pathways: methane metabolism, phenylalanine metabolism, and phenylpropanoid biosynthesis.
notably, the  <dig> genes enriched in the three pathways were all the same.
among them, four genes including glyma02g <dig>  glyma11g <dig>  glyma18g <dig> and glyma20g <dig> acted as peroxidase to generate p-hydroxyphenyl lignin , guaiacyl lignin  and syringyl lignin  from corresponding p-coumaryl alcohol, coniferyl alcohol and sinapyl alcohol .
one gene named glyma02g <dig> functioned in the production of p-coumaroyl-coa, which involved in the biosynthesis of flavonoid.fig.
<dig> common degs enriched on phenylpropanoid biosynthesis pathway in soybean
quantitative real-time pcr  verification
to verify the expression results obtained by gene chips, a total of  <dig> differently expressed genes either between soybean accessions cd and yh in low-p stress or between different p treatments in roots of cd and yh were selected for rt-qpcr.
the results of rt-qpcr analysis were consistent with the data obtained from gene chips , indicating the reliability of the results from gene chips.fig.
<dig> validation of the results obtained by gene chips with rt-qpcr
discussion
low phosphorus induced systemic changes in soybean
phosphorus is important to soybean growth and development.
to detect low-p-stress related genes, multiple comparative analyses were conducted.
first, we compared the expression level changes between different p treatments in the given materials and tissues, which represented intra-material and intra-tissue degs in coping with low-p conditions.
simultaneously, we analyzed inter-materials and inter-tissues gene expression patterns in the case of different concentrations of phosphorus.
finally, we obtained overlapping degs between inter-materials and intra-materials, between inter-tissues and intra-tissues, and common genes among different p treatments, materials and tissues.
through these comparative analyses, we acquired stable genes responding to low-p stress, and obtained common genes and pathways related to low-p stress in soybean.
zeng et al.
revealed  <dig> biological processes in go enrichment analysis of phosphate-deficiency-responsive genes in soybean roots   <cit> .
among these processes, photosynthesis, iron ion transport, fatty acid metabolic process, and stress responses were also detected in our study.
through the analyses of the go enrichment results, we found oxidative response, as well as peroxidases were repeatedly detected.
as earlier studies demonstrated, the expression levels of peroxidases were closely related to the resistance to biotic and abiotic stresses, such as salt stress, wounding, diseases, hypersensitive reactions, and so on .
recent studies confirmed that under p-limited conditions, alterations of oxidative stress- related genes including peroxidase genes were observed in zea mays, rice and arabidopsis  <cit> .
simultaneously, the results of go enrichment analyses also showed many genes induced by low-p stress were associated with plant hormone signal transduction, such as ethylene and jasmonic acid.
previous studies supported that many genes induced by low-p stress were associated with plant hormone signal transduction .
an increase of ethylene content in root systems of common bean indicated that ethylene regulated changes in root morphology  under low-p stress  <cit> .
meanwhile, the effect of ethylene on the elongation of primary roots in arabidopsis and the emergence of root hairs under low-p stress was shown with the use of a method based on image processing  <cit> .
gao et al.
studied the effect of ethylene on physiological changes of soybean plants grown in p-deficient solution.
their results suggested that ethylene increased root vigor and acid phosphatase activity in the root system, as well as root:shoot ratio under low-p stress  <cit> .
in addition, jasmonic acid promoted plants’ freezing tolerance with the application of exogenous jasmonic acid and blocking-up biosynthesis of jasmonic acid  <cit> .
the three common pathways, methane metabolism, phenylalanine metabolism, and phenylpropanoid biosynthesis, were enriched from comparative analyses, which suggested that these pathways responded actively when plants suffered low-p stress.
these results were further supported by the facts that transcript profiling analyses of zea mays and arabidopsis identified that many genes involved in biosynthesis of phenylpropanoids either up-regulated or down-regulated when exposed to low-p stress  <cit> .
phenylpropanoid was a member of phenolic compound, which pertained to plant secondary metabolites.
for example, hall et al.
found that the phenolic content was stimulated by an increase in nutrient stress in helianthus annuns  <cit> .
similar phenomenon occurred in tageteserecta when suffering from water stress with a significant increase of phenolic content  <cit> .
koeppe et al.
observed that more phenol compounds were leached from living intact roots, dried roots, and tops of phosphate-deficient plants than from phosphate-sufficient ones  <cit> .
although the effect of methane metabolism on the response to low-p stress had not been clearly stated, a recent study suggested plants released more methane at high temperatures  <cit> , and methane could play an important role in regulating plant growth  <cit> .
methane metabolism in association with low-p stress may provide a new direction in studying the mechanism of low-p tolerance in soybean.
our results suggested that the responses of soybean to low-p stress were complex cross-talks; it not only depended on other nutrition  and plant hormones levels , but also on the metabolism of secondary metabolites.
this information can be used to focus on these systemic connections induced by low-p stress to devise strategies aimed at improving soybean yield in p-deficiency soils.
massive genes involved in response to low-p stress in low-p-tolerant accession
we selected a low-p-tolerant accession, cd, and a low-p-sensitive accession, yh.
the numbers of different expression genes between different p treatments in cd were 30-fold more than those in yh.
to explore whether the degs were accession-specific between cd and yh, we compared the degs between these two accessions.
our results showed that among the  <dig> degs in yh, one gene  was also found in cd.
a homologous gene was identified in lotus japonicus  with the tool of blast in the national center for biotechnology information  database , which exhibited 54 % sequence identify to glyma20g <dig>  this gene was predicted as a nodulation-associated bzip transcription factor gene  <cit> .
for a deg glyma09g <dig> in yh, a homologous gene glyma17g <dig> was detected in the degs from cd.
glyma09g <dig> encoded a putative e <dig> ligase protein with metal ion binding and dna binding domains, homologous with brutus and oshrz in arabidopsis thaliana and oryza sativa which regulates the response to iron deficiency  <cit> .
the rest annotated degs in yh, three genes  encoded the transcription factor org2-like proteins.
the homologous gene in arabidopsis thaliana was bhlh <dig>  which played an important role in the iron-deficiency responses and uptake of this nutrient  <cit> .
glyma04g <dig> was a plastid-lipid-associated fibrillin protein, transferring phosphorus-containing groups.
glyma20g <dig> encoded a member of basic leucine zipper transcription gene family, homologs with tga <dig> in arabidopsis thaliana.
tga <dig> acts as a regulator in response to defense signals  <cit> .
although no homologous genes of the rest of the degs were detected in cd; gene response to iron-deficiency, uptake, and defense were detected in cd.
these genes might be located in different pathways or different locations of the same pathway.
in conclusion, most degs were specific between cd and yh , however, the involved pathways found in yh were included in cd.
in addition, other low-p stress related pathways were detected in cd, for example, oxidation-reduction and metabolic processes.
the much more genes induced by low-p stress existed in cd than in yh suggested massive genes involved in the response to low-p condition in tolerant accession cd.
in addition to involving in biotic and abiotic stress responses, some of these genes were related to plant growth and development.
on one hand, these massive genes could response to low-p condition through mediating the secretion of phytohormone and secondary metabolites.
on the other hand, they could improve uptake of p through controlling soybean root traits, such as main root elongation and increased root hair number.
all of these responses could make cd more adaptive to low-p stress and more resistant than yh.
likewise, more stress-responsive genes were identified in tolerant soybean line than those in sensitive soybean line in the response to common cutworm feeding  <cit> .
moreover, more degs were detected in cd than those in yh under low-p stress  could be due to different genetic backgrounds between these two materials, which was further supported by the result that there were  <dig>  snps detected between cd and yh by a genome-wide njau 355 k soysnp array .
the selection of p-efficient soybean accessions could be helpful in exploring low-p related genes and creating new low-p-tolerant soybean accessions.
roots played a crucial role in p metabolism during seedling stage
compared with  <dig> degs detected in roots under low-p stress, only  <dig> degs were detected in leaves .
plants under low-p conditions depended on the root system during seedling stage to absorb phosphorus nutrition.
previous studies showed an increase in lateral roots, root hair density, root:shoot ratio when plants were exposed to low-p stress  <cit> .
among them, uptake through root hairs contributed up to 90 % of the phosphorus acquired by plants  <cit> .
some phosphate transporters specifically located in roots had been identified to be involved in pi uptake at the root-soil interface in p-deficient soils.
these proteins transport available phosphorus that could be absorbed by plants, mainly h2po4−, hpo42−and po43−  <cit> .
the apparent changes of root systems under low-p stress in our hydroponic experiments also confirmed that roots played an important role in acquiring p during seedling stage.
as such, we found some genes associated with phosphorus absorption in roots.
the gene glyma04g <dig> showed homology with atspx2- <dig>  members of a sub-family of arabidopsis genes with the spx domain.
proteins harboring spx domain are believed to be involved in phosphorus acquisition and phosphorus signaling network.
duan et al.
found several pi starvation-responsive genes were regulated by atspx2- <dig>  positively or negatively.
further studies demonstrated the expression levels of atspx <dig> and atspx <dig> increased under phosphorus starvation  <cit> .
interestingly, we also found three degs in roots; namely glyma13g <dig> , glyma07g <dig>  and glyma13g <dig>  respectively, showed homology with pho <dig>  the ubiquitin-conjugating enzyme that contained spx domain, which is crucial for phosphorus acquisition and translocation in plants.
results of huang et al.
suggested pho <dig> modulated pi absorption, specifically transport at the root surface by regulating the abundance of pht1s  <cit> .
the limited number of degs in leaves probably resulted from genes in the leaves primarily involved in p distribution and utilization in later stages.
three degs from leaves  were homologous to bhlh <dig> in arabidopsis thaliana, functioning in the iron-deficiency responses and uptake  <cit> .
thus, we concluded roots played an important role in coping with low-p stress during seedling stage.
an in-depth study on degs in roots will advance the functional study of highly p-efficient genes.
CONCLUSIONS
p starvation leads to systemic changes in the gene expression of soybean.
here, the gene expression patterns of soybean under low-p stress were surveyed on the roots and leaves of a low-p-tolerant accession and a low-p-sensitive accession by microarray chips.
through comparative analyses on the differently expressed genes, we identified  <dig> candidate genes and three common pathways, including methane metabolism, phenylalanine metabolism and phenylpropanoid biosynthesis, which were highly correlated to low-p stress.
these results not only promote our understandings of the molecular bases of the responses to p deficiency, but also facilitate research in improving pi usage in soybean and designing highly phosphate-efficient soybeans.
this process could optimize fertilizer use and promote development of sustainable agricultural practices.

 as an essential constituent of many important compounds that are required for the development and growth of plants, phosphorus actively participates in photosynthesis, respiration, carbohydrate metabolism, energy transduction and other processes. in many agricultural systems, however, the concentration of absorbable phosphorus in soil is insufficient to ensure plant productivity. as one of the most important grain and industrial crops, soybean is a vital source of protein for human and animal food. however, soybean growth and yield are limited by varieties of edaphic factors, especially low phosphorus availability in soils. under low-p stress, soybean plants become dwarfed, leaf area is reduced, and necrotic spots might appear on lower leaves. also, flowering time is increased and the number of pods could decrease during the fruiting period. due to these p-mediated phenotypic effects on soybeans, phosphate fertilizers are used to overcome undesirable potential yield outcomes. large quantities of phosphate fertilizer may cause environmental pollution as it is easily immobilized in the soil to form chelates or precipitates. in order to alleviate the effects of phosphorus deficiency on soybean growth, it is particularly important for us to search for p-efficient soybean materials and investigate high-phosphorus-efficient genes. recently, four major groups of genes induced by low-p stress were reported. one group was involved in sensing low-p signals, including local phosphate sensing that was proposed to occur in the primary root in arabidopsis thaliana, and long-distance phosphate signaling which depends on the activity of transcription factors, such as phr phl. members of another group mediated phosphate distribution by involving an internal regulation of phosphorus homeostasis, including spx family and pho that contains spx domain in arabidopsis. the remaining two groups of genes mainly focus on the acquisition and transportation of phosphate. for example, gmacp encodes an acid phosphatase, the overexpression of which in roots could increase the absorption of phosphate in the nutrient solution. the pht family consisting of pht genes in soybean encoded phosphate transporters was demonstrated to be involved in not only directing phosphate uptake in roots, but also transporting phosphorus as phosphate. generally, the acquisition and transportation of phosphate greatly depends on the root system as root morphology changes under p-deficient conditions to sustain plant growth. this includes an increase in number and length of lateral roots in maize, longer and greater number of root hairs in low-p tolerant accessions of barley and arabidopsis, and relatively longer primary roots in rice when a low-phosphorus-related transcription factor osmyb2p- is over expressed. meanwhile, qtl mapping associated with soybean phosphorus efficiency have made processes. qtls related to seed phosphorus content, phosphorus efficiency relevant root traits, and soybean tolerance to low phosphorus stress based on flower and pod abscission rate have been identified. furthermore, screening of soybean p-efficient materials identified p-efficient genotypes, such as hp hp huaxai huaxia etc.. as an abiotic stress, low-p conditions may induce a series of complex responses. it was demonstrated that phytohormone signaling pathways responded strongly to phosphate deprivation. the expression of a large set of genes associated with plant hormone signaling either induced or repressed root growth in arabidopsis and zea mays when suffering low-p conditions. in addition, changes of the secondary metabolites to low-p stress were shown with an accumulation of anthocyanin under long-term phosphorus deficiency. similar results were observed in zea mays with genes related to the phenylpropanoid pathway identified as up-regulated or down-regulated under p starvation. a large set of genes and signaling pathways that were associated with low-p stress had been uncovered using the technique of gene expression profiles and transcriptome analyses in arabidopsis, rice, and zea mays. recently, zeng et al. identified phosphate-deficiency-responsive genes in soybean roots by high-throughput sequencing. selection of representative materials and consideration of different tissues will further reveal genes and pathways associated with low-p stress. here, we selected a low-p-tolerant accession chundou and a low-p-sensitive accession yunhefengwodou from soybean accessions through hydroponic experiments. then microarray chips were performed on the roots and leaves of these two accessions. by analyzing the gene expression data, we identified candidate genes and three common pathways that were induced by low-p stress. our study provides candidate loci for functional identification of high-phosphorus-efficient genes which may be of great significance for cultivating p-efficient soybean accessions.  identification of low-p-tolerant accession cd and low-p-sensitive accession yh with hydroponic experiments lasting for 10 days under different p treated conditions of soybean accessions, we found that the soybean accession cd and yh showed similar root morphology under + p treated condition, while showing obvious morphology differences under -p treated condition. we then identified accession cd as a low-p-tolerant accession, and accession yh as a low-p-sensitive accession. then, we carried out a gradient experiment using full, half and quarter hoagland’s solution and found that the + p and -p treated soybean plants showed most obvious morphology difference when grown in half hoagland. hence half hoagland with  mmol/l p and  mmol/l p were considered as normal and low p treated conditions, respectively. the phenotypic traits of cd and yh are shown in additional file 1: figure s primary root length was longer in cd than that in yh. the number of root hairs and lateral roots of cd were much more than those in yh. microarray experiments were performed on the roots and leaves of cd and yh. significantly differentially expressed genes in different p treated conditions under low p condition, and up-regulated differently expressed genes were found in the roots and leaves of low-p-tolerant accession cd respectively; while no up-regulated degs were found in the roots and leaves of low-p-sensitive accession yh compared with low-p treatment. for down-regulated degs, only degs could be found in the roots of cd, degs in the roots, and degs in the leaves of yh.table the number of degs in different p treated conditions generally, the number of degs in cd was more than that in yh, not only in the roots but also in the leaves, which suggested that cd could regulate a series of genes to cope with the low-p stress. we found only one common deg, glyma20g, in both the roots of cd and yh, with the fold change and down-regulated in low-p-tolerant accession cd and low-p-sensitive accession yh, respectively. the number of degs in the roots was more than that in the leaves. interestingly, one common deg, glyma01g named as gmmmp which encodes a matrix metalloproteinase, was found both in the roots and leaves of cd. the fold change is and, respectively. previous studies demonstrated that it was associated with pathogenic infections in soybean. a total of non-redundant degs were found between treated and control conditions in the roots and leaves of both soybean accessions used in this experiment. the hierarchical clustering analysis classified the tolerant and sensitive soybean accessions into two clusters; samples cluster and genes cluster. the compactness of the clusters of the three replicates confirmed the high reliability of our experimental design.fig. hierarchical clustering of differentially expressed genes in different p-treated conditions. a a and a three biological replicates of control groups; a a and a three biological replicates of treatment groups. red and green indicate high and low expression levels, respectively to evaluate the potential functions of the degs between low-p stress and normal conditions in both soybean accessions, gene ontology categories were applied for low-phosphorus-induced genes. the three go terms were categorized as biological processes, molecular functions, and cellular components. we found that degs were significantly associated with response to wounding, response to jasmonic acid stimulus, oxidoreductase activity, and cell wall, whether they were up-regulated genes or down-regulated genes. furthermore, we identified pathways that were associated with low-p induced genes with the use of kegg. overall, pathways with p values less than were motivated under low-p stress, and most pathways were related to metabolic process. notably, the most significant pathways that were enriched included phenylpropanoid biosynthesis, methane metabolism and phenylalanine metabolism, all related to the metabolism of secondary metabolites. significantly differentially expressed genes between different soybean accessions accessions cd and yh showed strong differences in the ability to suffer low-p stress. there were much more degs in cd than in yh. to further reveal the degs related with low-p condition between cd and yh, we compared the degs between soybean accessions cd and yh.table the number of degs in cd/yh and roots/leaves low p ↑a normal p → c normal p ↓ d a significantly up-regulated expression pattern under low-p stress; b significantly up-regulated expression pattern under normal p condition; c the expression patterns were not significantly changed under normal p condition; d significantly down-regulated expression pattern under normal p condition; e the number in parentheses indicates the number of genes identified by our criteria degs induced by low-p conditions were classified as three categories according to their expression patterns in normal conditions. for up-regulated genes in roots between cd and yh in low-p stress, the three categories were as follows: 1) genes were up-regulated in cd/yh in normal conditions, which is consistent with the expression pattern in low-p stress; 2) genes were down-regulated in normal conditions, which showed the opposite expression patterns between normal and low-p conditions; and 3) genes showed no expression difference in normal condition, with the expression difference less than two fold between cd and yh. by the criterion that genes with two fold expression change between experimental group and control group were identified as degs, the expression level in cd was more than two fold higher than that in yh for the up-regulated genes. among the up-regulated genes, there were genes with no expression difference in normal conditions, the expression difference of which between cd and yh was less than two fold. to make the expression difference of cd/yh between low-p stress and normal conditions higher than two fold, we chose those genes with larger than 4-fold expression difference between cd and yh under low-p conditions as candidates. by this criterion, genes with larger than 4-fold expression difference between cd and yh under low-p condition from the genes were chosen as candidates. hence, we got and degs from and degs, respectively. combined with the genes with opposite expression patterns, they were considered as marked degs between cd and yh. therefore, marked degs were obtained. among these genes, one gene was up-regulated more than fold between cd/yh under low-p conditions in both roots and leaves. thus, there were non-redundant marked degs detected in total. among these marked degs, the majority were up-regulated in cd/yh under low-p stress. additionally, there were and genes expressed differently between cd and yh in roots and leaves, respectively. this result suggested more genes were activated in roots than in leaves under low-p stress during seedling stage. we performed the enrichment of go and kegg analysis for these marked degs. go analysis showed that these genes were mainly involved in secondary metabolic processes, for example flavonoid biosynthesis. they were also involved in response to hormone stimuli, and defense and abiotic stimuli. in addition, they were involved in electron carrier activity and peroxidase. kegg analysis revealed the pathways enriched were mostly associated with the metabolism of secondary metabolites, including phenylpropanoid biosynthesis, methane metabolism, phenylalanine metabolism, flavonoid biosynthesis and so on. the degs in different p-treated conditions were considered as genes in response to low-p stress. the marked degs in different soybean accessions were not only induced by low-p stress but also by materials-variation dependent upon the response to low-p conditions. thus we compared degs and marked genes, finding overlapping genes that were considered as active genes induced by low-p stress. further study revealed most of the common genes focused on roots, which was consistent with a greater number of degs in roots than leaves in different p conditions. interestingly, the genes that showed opposite expression patterns in roots all belonged to the overlapping genes. despite two genes that could not find corresponding symbol and one gene that could not find functional annotation, seven genes with functional annotations were shown to be associated with biotic and abiotic stresses. likewise, go and kegg enrichment analyses were applied to evaluate the potential functions of the overlapping degs. degs significantly enriched in go terms of biological processes such as oxidation reduction and response to jasmonic acid stimulus. molecular functions such as oxidoreductase activity, and cellular components such as the cell wall, were also enriched. meanwhile, sixteen pathways were enriched based on significance level. similar to the results of the marked genes, the pathways of phenylpropanoid biosynthesis, methane metabolism, phenylalanine metabolism, and flavonoid biosynthesis were also detected. significantly differentially expressed genes between roots and leaves in our study, under low-p treatment, more degs were detected in roots than leaves, suggesting there were more genes in response to low-p stress during seedling stage in roots. degs between roots and leaves were found. degs induced by low-p conditions were classified into three categories according to their expression patterns in normal conditions. for up-regulated genes between roots and leaves in cd under low-p stress, three categories of degs were up-regulated genes, one down-regulated gene, and genes showing no different expression in roots/leaves in normal conditions. among the genes, there were genes that showed larger than 4-fold expression difference between roots and leaves under low-p conditions. additionally, we found and genes showing larger than 4-fold expression differences under low-p conditions from the and genes, respectively. with one opposite-expressed gene, there were marked degs between roots and leaves. we performed the enrichment of go and kegg analysis for these marked degs. go analysis showed these genes enriched in a group of biological processes, for example, oxidation reduction, response to oxidative stress, and response to abiotic stimulus. kegg analysis revealed most of these genes enriched in methane metabolism, phenylpropanoid biosynthesis, and phenylalanine metabolism. the marked degs between roots and leaves were related with low-p stress. there were degs in roots and leaves considered to be caused by low-p stress. thus, we compared degs and marked genes and found overlapping genes. interestingly, overlapping genes were up-regulated in cd between roots and leaves, consistent with the results that more genes were induced in cd than yh under low-p stress. go enrichment analysis showed the overlapping genes enriched in biological processes such as response to cadmium ion, defense, jasmonic acid stimulus. in the molecular function term and cellular component, these genes were significantly enriched in transcription regulator activity and in the nucleus. furthermore, we identified pathways including methane metabolism and flavonoid biosynthesis. common degs among treatments, materials, and tissues the marked degs between cd and yh under low-p stress were identified as degs between materials which were associated with low-p stress, while the marked degs between roots and leaves in low-p conditions were regarded as degs between tissues that were related to low-p stress. we compared the marked degs with the degs between different p treatments and common degs were found to be in all three classes. that is to say, the degs occurred not only between treatments, but also between materials and tissues under low-p stress, which suggested they were active genes under low-p conditions. further kegg analyses enriched in three general pathways: methane metabolism, phenylalanine metabolism, and phenylpropanoid biosynthesis. notably, the genes enriched in the three pathways were all the same. among them, four genes including glyma02g glyma11g glyma18g and glyma20g acted as peroxidase to generate p-hydroxyphenyl lignin, guaiacyl lignin and syringyl lignin from corresponding p-coumaryl alcohol, coniferyl alcohol and sinapyl alcohol. one gene named glyma02g functioned in the production of p-coumaroyl-coa, which involved in the biosynthesis of flavonoid.fig. common degs enriched on phenylpropanoid biosynthesis pathway in soybean quantitative real-time pcr verification to verify the expression results obtained by gene chips, a total of differently expressed genes either between soybean accessions cd and yh in low-p stress or between different p treatments in roots of cd and yh were selected for rt-qpcr. the results of rt-qpcr analysis were consistent with the data obtained from gene chips, indicating the reliability of the results from gene chips.fig. validation of the results obtained by gene chips with rt-qpcr discussion low phosphorus induced systemic changes in soybean phosphorus is important to soybean growth and development. to detect low-p-stress related genes, multiple comparative analyses were conducted. first, we compared the expression level changes between different p treatments in the given materials and tissues, which represented intra-material and intra-tissue degs in coping with low-p conditions. simultaneously, we analyzed inter-materials and inter-tissues gene expression patterns in the case of different concentrations of phosphorus. finally, we obtained overlapping degs between inter-materials and intra-materials, between inter-tissues and intra-tissues, and common genes among different p treatments, materials and tissues. through these comparative analyses, we acquired stable genes responding to low-p stress, and obtained common genes and pathways related to low-p stress in soybean. zeng et al. revealed biological processes in go enrichment analysis of phosphate-deficiency-responsive genes in soybean roots. among these processes, photosynthesis, iron ion transport, fatty acid metabolic process, and stress responses were also detected in our study. through the analyses of the go enrichment results, we found oxidative response, as well as peroxidases were repeatedly detected. as earlier studies demonstrated, the expression levels of peroxidases were closely related to the resistance to biotic and abiotic stresses, such as salt stress, wounding, diseases, hypersensitive reactions, and so on. recent studies confirmed that under p-limited conditions, alterations of oxidative stress- related genes including peroxidase genes were observed in zea mays, rice and arabidopsis. simultaneously, the results of go enrichment analyses also showed many genes induced by low-p stress were associated with plant hormone signal transduction, such as ethylene and jasmonic acid. previous studies supported that many genes induced by low-p stress were associated with plant hormone signal transduction. an increase of ethylene content in root systems of common bean indicated that ethylene regulated changes in root morphology under low-p stress. meanwhile, the effect of ethylene on the elongation of primary roots in arabidopsis and the emergence of root hairs under low-p stress was shown with the use of a method based on image processing. gao et al. studied the effect of ethylene on physiological changes of soybean plants grown in p-deficient solution. their results suggested that ethylene increased root vigor and acid phosphatase activity in the root system, as well as root:shoot ratio under low-p stress. in addition, jasmonic acid promoted plants’ freezing tolerance with the application of exogenous jasmonic acid and blocking-up biosynthesis of jasmonic acid. the three common pathways, methane metabolism, phenylalanine metabolism, and phenylpropanoid biosynthesis, were enriched from comparative analyses, which suggested that these pathways responded actively when plants suffered low-p stress. these results were further supported by the facts that transcript profiling analyses of zea mays and arabidopsis identified that many genes involved in biosynthesis of phenylpropanoids either up-regulated or down-regulated when exposed to low-p stress. phenylpropanoid was a member of phenolic compound, which pertained to plant secondary metabolites. for example, hall et al. found that the phenolic content was stimulated by an increase in nutrient stress in helianthus annuns. similar phenomenon occurred in tageteserecta when suffering from water stress with a significant increase of phenolic content. koeppe et al. observed that more phenol compounds were leached from living intact roots, dried roots, and tops of phosphate-deficient plants than from phosphate-sufficient ones. although the effect of methane metabolism on the response to low-p stress had not been clearly stated, a recent study suggested plants released more methane at high temperatures, and methane could play an important role in regulating plant growth. methane metabolism in association with low-p stress may provide a new direction in studying the mechanism of low-p tolerance in soybean. our results suggested that the responses of soybean to low-p stress were complex cross-talks; it not only depended on other nutrition and plant hormones levels, but also on the metabolism of secondary metabolites. this information can be used to focus on these systemic connections induced by low-p stress to devise strategies aimed at improving soybean yield in p-deficiency soils. massive genes involved in response to low-p stress in low-p-tolerant accession we selected a low-p-tolerant accession, cd, and a low-p-sensitive accession, yh. the numbers of different expression genes between different p treatments in cd were 30-fold more than those in yh. to explore whether the degs were accession-specific between cd and yh, we compared the degs between these two accessions. our results showed that among the degs in yh, one gene was also found in cd. a homologous gene was identified in lotus japonicus with the tool of blast in the national center for biotechnology information database, which exhibited 54 % sequence identify to glyma20g this gene was predicted as a nodulation-associated bzip transcription factor gene. for a deg glyma09g in yh, a homologous gene glyma17g was detected in the degs from cd. glyma09g encoded a putative e ligase protein with metal ion binding and dna binding domains, homologous with brutus and oshrz in arabidopsis thaliana and oryza sativa which regulates the response to iron deficiency. the rest annotated degs in yh, three genes encoded the transcription factor org2-like proteins. the homologous gene in arabidopsis thaliana was bhlh which played an important role in the iron-deficiency responses and uptake of this nutrient. glyma04g was a plastid-lipid-associated fibrillin protein, transferring phosphorus-containing groups. glyma20g encoded a member of basic leucine zipper transcription gene family, homologs with tga in arabidopsis thaliana. tga acts as a regulator in response to defense signals. although no homologous genes of the rest of the degs were detected in cd; gene response to iron-deficiency, uptake, and defense were detected in cd. these genes might be located in different pathways or different locations of the same pathway. in conclusion, most degs were specific between cd and yh, however, the involved pathways found in yh were included in cd. in addition, other low-p stress related pathways were detected in cd, for example, oxidation-reduction and metabolic processes. the much more genes induced by low-p stress existed in cd than in yh suggested massive genes involved in the response to low-p condition in tolerant accession cd. in addition to involving in biotic and abiotic stress responses, some of these genes were related to plant growth and development. on one hand, these massive genes could response to low-p condition through mediating the secretion of phytohormone and secondary metabolites. on the other hand, they could improve uptake of p through controlling soybean root traits, such as main root elongation and increased root hair number. all of these responses could make cd more adaptive to low-p stress and more resistant than yh. likewise, more stress-responsive genes were identified in tolerant soybean line than those in sensitive soybean line in the response to common cutworm feeding. moreover, more degs were detected in cd than those in yh under low-p stress could be due to different genetic backgrounds between these two materials, which was further supported by the result that there were snps detected between cd and yh by a genome-wide njau 355 k soysnp array. the selection of p-efficient soybean accessions could be helpful in exploring low-p related genes and creating new low-p-tolerant soybean accessions. roots played a crucial role in p metabolism during seedling stage compared with degs detected in roots under low-p stress, only degs were detected in leaves. plants under low-p conditions depended on the root system during seedling stage to absorb phosphorus nutrition. previous studies showed an increase in lateral roots, root hair density, root:shoot ratio when plants were exposed to low-p stress. among them, uptake through root hairs contributed up to 90 % of the phosphorus acquired by plants. some phosphate transporters specifically located in roots had been identified to be involved in pi uptake at the root-soil interface in p-deficient soils. these proteins transport available phosphorus that could be absorbed by plants, mainly h2po4−, hpo42−and po43−. the apparent changes of root systems under low-p stress in our hydroponic experiments also confirmed that roots played an important role in acquiring p during seedling stage. as such, we found some genes associated with phosphorus absorption in roots. the gene glyma04g showed homology with atspx2- members of a sub-family of arabidopsis genes with the spx domain. proteins harboring spx domain are believed to be involved in phosphorus acquisition and phosphorus signaling network. duan et al. found several pi starvation-responsive genes were regulated by atspx2- positively or negatively. further studies demonstrated the expression levels of atspx and atspx increased under phosphorus starvation. interestingly, we also found three degs in roots; namely glyma13g, glyma07g and glyma13g respectively, showed homology with pho the ubiquitin-conjugating enzyme that contained spx domain, which is crucial for phosphorus acquisition and translocation in plants. results of huang et al. suggested pho modulated pi absorption, specifically transport at the root surface by regulating the abundance of pht1s. the limited number of degs in leaves probably resulted from genes in the leaves primarily involved in p distribution and utilization in later stages. three degs from leaves were homologous to bhlh in arabidopsis thaliana, functioning in the iron-deficiency responses and uptake. thus, we concluded roots played an important role in coping with low-p stress during seedling stage. an in-depth study on degs in roots will advance the functional study of highly p-efficient genes.  p starvation leads to systemic changes in the gene expression of soybean. here, the gene expression patterns of soybean under low-p stress were surveyed on the roots and leaves of a low-p-tolerant accession and a low-p-sensitive accession by microarray chips. through comparative analyses on the differently expressed genes, we identified candidate genes and three common pathways, including methane metabolism, phenylalanine metabolism and phenylpropanoid biosynthesis, which were highly correlated to low-p stress. these results not only promote our understandings of the molecular bases of the responses to p deficiency, but also facilitate research in improving pi usage in soybean and designing highly phosphate-efficient soybeans. this process could optimize fertilizer use and promote development of sustainable agricultural practices. 
the functional classification of these genes revealed gene ontology terms of biological process including response to oxidation reduction, hormone stimuli, and biotic and abiotic stimuli. genome-wide analyses were performed on the roots and leaves of a low-p-tolerant accession and a low-p-sensitive accession which were identified by hydroponic experiments under different p treatments. additionally, three common pathways were identified. plants have developed strategies to cope with this, while few genes associated with low-p tolerance have been identified in soybean. these results could not only promote the work on the molecular regulation mechanism under low-p stress in soybean, but also facilitate the cultivation of high-phosphorus-acquisition and high-phosphorus-utilization soybean varieties. through comparative analyses on the differently expressed genes, we explored common genes that were highly correlated to low-p stress. however, stress under low-p conditions has become a limiting factor that affects crop yields and qualities. soybeanroot systemlow-p stressphosphorus deficiencyexpression profilekey transgenic breeding program of china2014zx08004-003yu deyue http://dx.doi.org/ /501100001809national natural science foundation of china 31301342wang jiao http://dx.doi.org/ /501100001809national natural science foundation of china31370034cheng hao http://dx.doi.org/ /501100001809national natural science foundation of china 31301336zhang dan jiangsu collaborative innovation center for modern crop productionissue-copyright-statement© the author 2016
 phosphorus is one of the most important macronutrients that is required for plant growth and development. however, stress under low-p conditions has become a limiting factor that affects crop yields and qualities. plants have developed strategies to cope with this, while few genes associated with low-p tolerance have been identified in soybean.  genome-wide analyses were performed on the roots and leaves of a low-p-tolerant accession and a low-p-sensitive accession which were identified by hydroponic experiments under different p treatments. through comparative analyses on the differently expressed genes, we explored common genes that were highly correlated to low-p stress. the functional classification of these genes revealed gene ontology terms of biological process including response to oxidation reduction, hormone stimuli, and biotic and abiotic stimuli. additionally, three common pathways were identified.  these results could not only promote the work on the molecular regulation mechanism under low-p stress in soybean, but also facilitate the cultivation of high-phosphorus-acquisition and high-phosphorus-utilization soybean varieties. electronic supplementary material the online version of this article contains supplementary material, which is available to authorized users. keywords soybeanroot systemlow-p stressphosphorus deficiencyexpression profilekey transgenic breeding program of china2014zx08004-003yu deyue http://dx.doi.org/ /501100001809national natural science foundation of china 31301342wang jiao http://dx.doi.org/ /501100001809national natural science foundation of china31370034cheng hao http://dx.doi.org/ /501100001809national natural science foundation of china 31301336zhang dan jiangsu collaborative innovation center for modern crop productionissue-copyright-statement© the author 2016 