Genetic material
A population composed of 396 recombinant inbred lines (RILs, F10:11), derived from a cross between the sorghum lines, BR007 and SC283, was developed by single-seed descent [47] at Embrapa Maize and Sorghum (Sete Lagoas - MG, Brazil). Both BR007 (Redbine-type) and SC283 (sorghum converted guinea) were introduced into the Embrapa breeding program in 1972 from the Purdue Breeding Program (West Lafayette - IN, US). BR007 is Al sensitive whereas SC283 is highly tolerant to Al toxicity [33]. Previous studies indicated that, while SC283 has higher grain yield in a soil with low-P availability compared to BR007, the grain yield increase in BR007 in response to adequate P supply is in turn higher than in SC283 [48].
Phenotyping for low-P in field conditions
Four field experiments were conducted at the experimental station of Embrapa Maize and Sorghum in Sete Lagoas, State of Minas Gerais, Brazil, during the summer season of 2012–2013. The experimental site is a clay and highly weathered tropical soil, with low fertility in natural conditions, low pH, Al toxicity and low-P. Soil P (Mehlich 1) varied from 1 to 6 ppm between 0 and 20 cm of soil depth, and from 1 to 4 ppm at the sub-superficial soil layer (20–40 cm). The minimum and maximum content of available P in the soil (Psoil) was 5.88 kg ha− 1 and 19.79 kg ha− 1.
Each experiment was arranged as a 12 × 10 alpha lattice design, with three complete replicates and ten incomplete blocks per replicate. Each block contained 12 plots, within which ten RILs (regular treatments) and the two parents (common checks) were allocated. Each plot consisted of a three-meter row, with 0.45 m between rows and 8 plants m− 1. Fertilization consisted of 150 kg ha− 1 of 20–00-20 (NPK) at sowing and 200 kg ha− 1 of urea applied 30 days after sowing.
Grain yield (Gy, kg ha− 1), flowering time (FT, days), plant height (PH, cm), phosphorus content in the plant (leaves and stems - Pp, kg ha− 1) and phosphorus content in the grain (Pg, kg ha− 1) were evaluated. For P measurements, samples of plant tissues and grains were collected in each plot, weighted and then dried at 65 °C to constant weight. Dry plant tissues and grains were then weighted, grounded and homogenized. Twenty gram - subsamples were used to determine P concentration and total P content (Pt), using inductively-coupled argon plasma emission spectrometry.
The phosphorus efficiency indexes were calculated according to the methodology proposed by Moll et al. [32], where: 1) phosphorus use efficiency (PUE) is equal to the product between phosphorus acquisition efficiency (PAE) and phosphorus internal utilization efficiency (PUTIL); 2) PAE is the total phosphorus content (Pt = Pp, P content in the plant + Pg, P content in the grain) divided by P content available in the soil; 3) PUTIL is Gy divided by Pt.
Root system phenotyping in low-P conditions
Root morphology traits were assessed in nutrient solution as described by de Sousa et al. [49] and Hufnagel et al. [22], using a randomized block design with three replicates. Seeds were surface-sterilized using sodium hypochlorite (5%), washed with distilled water and placed in moistened paper rolls. After 4 days, uniform seedlings were transferred to moistened blotting papers and placed into paper pouches (24 × 33 × 0.02 cm) as described by Hund et al. [50].
Each experimental unit consisted of one pouch, with three plants per pouch, whose bottom (3 cm) was immersed in containers filled with 5 l of the nutrient solution described in [51], with pH 5.65 and a P concentration of 2.5 μM. The containers were kept in a growth chamber with 27 °C day and 20 °C night temperatures and a 12-h photoperiod, under continuous aeration for 13 days.
After 13 days, root images were acquired using a digital camera Nikon D300S SLR. Images were then analyzed using both the RootReader2D (http://www.plantmineralnutrition.net/software/rootreader2d/) and WinRhizo (http://www.regent.qc.ca/) software. The following traits were measured: root length (RL - cm); root diameter (RD - mm); total root surface area (SA - cm2); surface area of very fine roots between 0 and 1 mm in diameter (SA1 - cm2); surface area of fine roots between 1 and 2 mm in diameter (SA2 - cm2); surface area of thicker roots between 2 and 4.5 mm in diameter (SA3 - cm2); root volume (RV - cm3) and volume of fine roots between 1 and 2 mm in diameter (V2 - cm3). Shoot dry matter (SDM) and root dry matter (RDM), phosphorus content in the shoot (Ps) and phosphorus content in the root (Pr) (in grams) were also measured.
Phenotypic analyses
Traits assessed in field and hydroponic experiments were analyzed using mixed models. For field experiments, the following model was used:
$$ {y}_{ijkl}=\mu +{E}_j+{R}_{k(j)}+{B}_{l(kj)}+{G}_i+{\varepsilon}_{ijkl} $$
yijkl is the phenotypic value of individual i in the block l of the kth replicate, within the experiment j; μ is the overall mean; and Gi is the genetic effect of individual i, which can be defined as:
$$ {G}_i=\left\{\begin{array}{c}{g}_i\ i=1,\dots, {n}_g\ \\ {}{t}_ii={n}_g+1,\dots, {n}_g+{n}_c\end{array}\right. $$
gi is the random effect of RIL i, ng is the total number of RILs; ti is the fixed effect of check i; and nc is the total number of checks; Ej is the fixed effect of the jth experiment (j = 1, … , 4); Rk(j) is the fixed effect of replicate k (k = 1, … , 3) in experiment j; Bl(kj) is the random effect of block l (l = 1, … , 10) in the replicate k, within the experiment j; and ε = (ε1111, ε2111, … , εIJKL)′ is a Nobs × 1 residual random vector assumed to be normally distributed with mean zero and variance \( {\sigma}_{\varepsilon}^2 \), in which Nobs is the total number of observations.
The model used for analyzing the hydroponic experiments was:
$$ {y}_{ij}=\mu +{B}_j+{g}_i+{\varepsilon}_{ij} $$
where yij is the phenotypic value of the RIL i (i = 1, … , ng) in the block j; μ is the overall mean; gi is the random genetic effect of RIL i; Bj is the fixed effect of block j (j = 1, … , 3); and ε = (ε11, ε21, … , εIJ)′ is a Nobs × 1 residual random vector assumed to be normally distributed with mean zero and variance \( {\sigma}_{\varepsilon}^2 \). Fixed and random effects were tested using the Wald statistics [52] and the likelihood ratio test (LRT, [53]) respectively, considering a 5% significance level (α).
For both statistical models, the genetic effect of RIL was first taken as random for estimating the genetic variance component (\( {\sigma}_g^2 \)) via restricted maximum likelihood (REML), and the heritability coefficient of each trait. The effect of RIL was then considered as fixed for estimating the adjusted means using best linear unbiased estimators (BLUEs). All the mixed models analyses were performed using the GenStat software (v.17.1.0) [54].
Trait heritabilities were estimated as proposed by Cullis et al. [55], called generalized heritabilities, using:
$$ {h}^2=1-\frac{\overline{v} BLUP}{2{\sigma}_g^2} $$
where \( \overline{v} BLUP \) is the average variance of the difference between two best linear unbiased predictions (BLUPs). Person’s correlation coefficients [56] were estimated based on the adjusted means of genotypes for traits assessed in the field and in the hydroponic experiments, using the package Hmic [57] in R [58].
SNP markers
Genomic DNA was isolated from approximately 500 mg of leaf tissue (eight plants per accession, i.e. RILs and their parents) as described by Saghai-Maroof et al. [59]. DNA samples were genotyped by sequencing according to Elshire et al. [60]. Reads were aligned to the version 1.4 of sorghum reference genome using Burrows-Wheeler Aligner program (BWA - [61]), and the SNP calling was performed using the GBS pipeline [62] implemented in the TASSEL software [63]. Missing genotypes were imputed using the NPUTE software [64]. Then, SNP data were filtered for 40% of minor allele frequency (MAF).
QTL mapping
The final set of traits used for multi-trait QTL mapping was comprised of grain yield (Gy), surface area of fine roots in the 1–2 mm diameter class (SA2) and root diameter (RD). Multi-trait QTL mapping analysis was performed according to the procedures described in Silva et al. [65] implemented in R. For that, a multi-locus QTL mapping procedure was considered, using the Haley & Knott regression [66] and the following linear model:
$$ {y}_{ti}={\mu}_t+\sum \limits_{r=1}^m{a}_{tr}{x}_{ir}+{\varepsilon}_{ti} $$
where yti is the adjusted mean of RIL i (i = 1, … , ng) for trait t (t = 1, … , T); μt is the intercept for each trait; atr is the rth QTL main effect on trait t; xir represents the genotype of RIL i for the SNP marker r (r = 1, … , nM), being nM the total number of markers; xir assumed values equal to 0 or 2 for RILs with homozygous genotypes for the allele donated by BR007B or SC283, respectively; and εi = (ε1i, ε2i, … , εTi)′ is a T × 1 random vector assumed to be independent and identically distributed according to a multivariate normal distribution with mean vector zero and positive definite symmetric variance-covariance matrix Σε, i.e. εi~MVN(0, Σε). Single-trait QTL mapping analyses were performed for each trait, using the above model (t = 1 gives an univariate regression model).
Multiple QTL models were built based on a forward-selection procedure, testing the significance of a putative QTL main effect at each SNP position along the genome. Significance of QTLs main effects were tested using the Score Statistic [67], considering a 10% significance level (α). According to simulations performed by Silva et al. [65], this significance level maximized the QTL detection power and kept the false discovery rate (i.e. the proportion of spurious QTLs) within an acceptable level.
QTL positions were refined after the inclusion of every new QTL in the model, until no more significant QTL main effects were found. Finally, non-significant QTL effects were removed from the model in a backward-elimination procedure, such as proposed in the seemingly unrelated regression coefficients method [68], considering a 1% significance level.
Quantitative analysis of SbPSTOL1 gene expression
Sorghum seedlings were grown in a modified Magnavaca nutrient solution [51] containing a low-P concentration (2.5 μM P), as described in the section Root system phenotyping in low-P conditions. The experiment was set up in randomized block design with three replicates and three seedlings per experimental unit (paper pouch), giving a total of nine biological replicates per genotype. After 13 days in nutrient solution, the expression profiles of the SbPSTOL1-like genes (Sb03g006765, Sb03g031690, and Sb07g002840) were assessed in the roots of the RIL parents, BR007 and SC283. Total RNA was isolated from bulked root tissues (nine roots per bulk), using the SV Total RNA Isolation System kit (Promega Corporation, Madison, WI, USA), according to the manufacturer’s instructions. Total RNA (1 μg) was used for cDNA synthesis using the High Capacity cDNA Reverse Transcription kit (Applied Biosystems, Foster City, CA, USA). Transcripts were quantified by quantitative real-time PCR (qPCR-RT), using SYBR Green technology with the ABI Prism 7500 Fast System (Applied Biosystems, Foster City, CA, USA).
Transcript relative quantification was performed with 20 ng cDNA samples and 0.02 ng for the endogenous constitutive control (18 s rRNA). Primers were designed for SbPSTOL1 and 18 s rRNA sorghum genes using the PrimerQuest tool (https://www.idtdna.com/PrimerQuest/) (Additional file 7). Calculation of relative gene expressions were performed using the 2-ΔΔCT method [69], with three technical replicates.