Major Novel QTL for Resistance to Cassava Bacterial Blight Identified through a Multi-Environmental Analysis.

Cassava, Manihot esculenta Crantz, has been positioned as one of the most promising crops world-wide representing the staple security for more than one billion people mainly in poor countries. Cassava production is constantly threatened by several diseases, including cassava bacterial blight (CBB) caused by Xanthomonas axonopodis pv. manihotis (Xam), it is the most destructive disease causing heavy yield losses. Here, we report the detection and localization on the genetic map of cassava QTL (Quantitative Trait Loci) conferring resistance to CBB. An F1 mapping population of 117 full sibs was tested for resistance to two Xam strains (Xam318 and Xam681) at two locations in Colombia: La Vega, Cundinamarca and Arauca. The evaluation was conducted in rainy and dry seasons and additional tests were carried out under controlled greenhouse conditions. The phenotypic evaluation of the response to Xam revealed continuous variation. Based on composite interval mapping analysis, 5 strain-specific QTL for resistance to Xam explaining between 15.8 and 22.1% of phenotypic variance, were detected and localized on a high resolution SNP-based genetic map of cassava. Four of them show stability among the two evaluated seasons. Genotype by environment analysis detected three QTL by environment interactions and the broad sense heritability for Xam318 and Xam681 were 20 and 53%, respectively. DNA sequence analysis of the QTL intervals revealed 29 candidate defense-related genes (CDRGs), and two of them contain domains related to plant immunity proteins, such as NB-ARC-LRR and WRKY.

Introduction

Cassava, Manihot esculenta Crantz, is a starchy root crop and one of the main staple food crops over the world due to its essential role for food security in tropical regions. This crop represents an important source of calories for about one billion people (Ceballos et al., 2010). Cassava tolerates drought, therefore it has been considered as one of the best alternatives for providing food for the world population in the context of climatic change (Howeler et al., 2013). The major bacterial vascular disease affecting this crop is Cassava Bacterial Blight (CBB), caused by Xanthomonas axonopodis pv. manihotis (Xam). The disease has a very high destructive power causing losses between 12 and 100% in affected areas (Lozano, 1986; López and Bernal, 2012). Xam has been described among the top 10 most important plant pathogenic bacteria (Mansfield et al., 2012). CBB has been reported in all regions where cassava is grown (López and Bernal, 2012; Taylor et al., 2012), including 56 countries distributed over Asia, Africa, Oceania and North, Central, and South America (http://www.cabi.org/) and the number of countries affected by the disease is increasing. The advancement of CBB epidemics has been reported in many countries, with Burkina Faso, being one of the most recent ones (Wonni et al., 2015). The Colombian Xam populations remain highly dynamic and exhibit a high genetic diversity (Trujillo et al., 2014). The analysis of 65 Xam genomes revealed that this pathogen harbors 14–22 effector genes, from which nine are conserved in all the strains (Bart et al., 2012).

Although some measures such as planting of disease-free material can be applied in the cassava fields in order to protect the crop, the safest and most efficient strategy to control CBB is to take advantage of natural plant genetic resistance to develop resistant cultivars for cultivation in CBB prone regions. Plants have evolved several mechanisms to defend themselves against pathogens. These mechanisms have been extensively studied in model plants, but knowledge generated in cassava is relatively scarce. Histology and cytochemistry studies of the resistance mechanisms in cassava during Xam infection showed callose deposition that act as a barrier in cortical parenchyma cells and phloem to block bacterial multiplication and dispersion (Kpémoua et al., 1996; Sandino et al., 2015). Other mechanisms of defense response including cell wall fortification, lignification and suberization associated with callose deposition and production of flavonoids and polysaccharides were also observed and they are faster and stronger in resistant cultivars compared to susceptible ones (Kpémoua et al., 1996). On the other side, efforts have been made in the last years to identify molecular determinants of the CBB resistance, including the selective analysis of homologous genes coding for proteins containing NBS and TIR domains (Lopez et al., 2003) or by annotation of sequence information based on the cassava genome sequence draft (Lozano et al., 2015; Soto et al., 2015).

Resistance to CBB has been described as quantitative, showing polygenic and additive inheritance (Hahn et al., 1974; Jorge et al., 2000, 2001) occurring together with resistance to cassava mosaic disease (CMD) (Hahn et al., 1980; Lokko et al., 2004; Rabbi et al., 2014). Several quantitative trait loci for resistance to CBB have been identified in cassava using full-sib population derived from the cross TMS30572 × CM1477-2. Eight QTL, explaining between 7.2 and 18.2% of the variance were detected in field conditions under high disease pressure and over two consecutive crop cycles (Jorge et al., 2001). Under controlled conditions, 12 resistance QTL were identified to five Xam strains, explaining 9–27% of the phenotypic variance (Jorge et al., 2000). Two other QTL were identified for resistance against Xam strains CIO151 and CIO121 explaining 62 and 21% resistance, respectively (Lopez et al., 2007). Moreover, Wydra et al. (2004) reported nine QTL explaining from 16 to 33% of the resistance variance to four African Xam strains.

The environment plays an important role in the phenotypic variation of quantitative traits (Weinig and Schmitt, 2004; Anderson et al., 2014). Consequently, the detection and stability of the QTL between environments is an important aspect to consider when studying genetic determinants of complex traits (Anderson et al., 2013; Mitchell-Olds, 2013; El-Soda et al., 2014). A QTL by environment interaction (Q × E) is defined as a different expression of the trait in varying environments, like having a significant effect in one environment but not in another (El-Soda et al., 2014). The complexity of the analysis increases when the plant-pathogen interaction under various environmental conditions is studied. According to the classic quantitative genetic model, the phenotype is the result of the genetic composition of the plant (G), the environment (E) and of the interaction between them (GxE). However, genotypes by environment interactions are additionally complicated in the case of plant-pathogen interactions in which another genotype (corresponding to the pathogen) must be considered. In this case, the equation of this form of interaction will result in interaction = Gplant × Gpathogen × E (Jorgensen, 2012). In this context, the detection of QTL and the dissection of their allelic composition are necessary to better understand the genetic basis of these interactions and the phenotypic responses to specific environmental conditions (El-Soda et al., 2014).

The molecular bases and mechanisms of quantitative resistance in plants are not well known. Interestingly, recent reports, including the cloning of genes operating in QTL, have provided some important clues that would facilitate the unraveling of the molecular bases of quantitative resistance in cassava. The proteins involved in pathogen recognition can be similar in both qualitative and quantitative resistance (Gebhardt and Valkonen, 2001; Lopez et al., 2003; Ramalingam et al., 2003; Calenge and Durel, 2006; Poland et al., 2009; Kou and Wang, 2010; Roux et al., 2014; Corwin et al., 2016). In addition, genes involved in quantitative resistance can correspond to those coding for proteins involved in the signal transduction pathways of defense response (Fukuoka et al., 2009; Corwin et al., 2016) and/or possessing antimicrobial activity (Ramalingam et al., 2003; Guimaraes and Stotz, 2004; Liu et al., 2006; Van Loon et al., 2006).

The aim of this study was to better understand the mechanistic basis of plant response to CBB. By analyzing a segregating population for disease resistance in several locations with various environmental conditions allows not only the identification of genetic factors leading to resistance, but as well the interaction of specific genotype with the environment.

Here, we report the identification of 5 novel QTL for specific resistance against two Xam strains in cassava, detected in field evaluations during rainy and dry seasons in two Colombian environments, as well under greenhouse conditions. Furthermore, the analysis of G × E and QTL × E interactions (Q × E), allowed to estimate the impact of the environment over the QTL effect. Candidate defense-related genes (CDRGs) located within the QTL intervals were identified based on the in silico analysis of DNA sequence information from the corresponding chromosome regions using whole genome sequence draft of cassava.

Materials and methods

Plant material

The mapping population consists of 117 full sib F1 segregating population derived from a cross between Nigerian cultivar TMS30572 and CIAT's elite cultivar CM2177-2 (Fregene et al., 1997). This population has been used extensively in mapping studies (Fregene et al., 1997; Jorge et al., 2000; Mba et al., 2001; Lopez et al., 2007; Soto et al., 2015). To produce multiple clonal plants of each individual, parent and F1 progeny were grown and propagated vegetatively from stem-cuttings in CBB-free field conditions at La Vega (Cundinamarca) and Universidad Nacional de Colombia, at Arauca (Arauca).

Bacterial strains and phenotypic evaluation of CBB resistance

For long-term storage, bacteria strains Xam318 and Xam681 were kept in 60% glycerol at −80°C, and streaked on LPGA medium [yeast extract (5 g/Lt), peptone (5 g/Lt), glucose (5 g/Lt), bacto-agar (15 g/Lt)] at 28°C for 12 h before use as inoculum. A single colony of Xam was grown in liquid culture at 28°C with shaking at 230 rpm for 24 h. Cells were harvested by centrifugation at 3,000 g and re-suspended in 10 mM MgCl2. Based on a preliminary screen, the Xam strains Xam318 and Xam681 were selected because the parent TMS30572 is resistant to these two strains while line CM2177-2, is susceptible (data not shown and Trujillo et al., 2014). The Xam strains were isolated from plant samples with typical symptoms of CBB disease in fields located at Ciénaga de Oro (N 08.889°, W 075.569°) a typical Caribbean savanna region and Palmitos (N 09.450°, W 075,160°) located on mountains.

For the evaluation of CBB disease response, plants were grown on two different locations: La Vega Cundinamarca, Andina region (Latitude 0,5°,00′44.188″N, Longitude 74°21′31.005″N) and Arauca, Arauca, Orinoquía region (Latitude 7° 1′ 22.32″N, Longitude 70° 44′ 42.50″W). Six-week old plants were inoculated on July 2013 and December 2014, corresponding to rainy and dry seasons, respectively. Besides, an evaluation was done under greenhouse controlled conditions. The corresponding codes for the different conditions, environments, strains and seasons, are presented in Table 1.

Table 1

Location, season, temperatures, relative humidity and mean precipitation of the environments, Xam strain and code where the inoculation and phenotyping were conducted.

Location	Season	Max. Temperature (C°)a	Min. Temperature (C°)a	Relative humidity (%)a	Mean precipitation (mm)a	Xam strain	Code
Arauca (Arauca)	Rainy	31	22	88	301	Xam318	AR318-R
	Dry	31	22	73	18,7	Xam318	AR318-D
	Rainy	31	22	88	301	Xam681	AR681-R
	Dry	31	22	73	18,7	Xam681	AR681-D
La Vega (Cundinamarca)	Rainy	20	9	79	106	Xam318	LV318-R
	Dry	22	10	70	30	Xam318	LV318-D
	Rainy	20	9	79	106	Xam681	LV681-R
	Dry	22	10	70	30	Xam681	LV681-D
Greenhouse 12 h of photoperiod	2013	30	20	70	–	Xam318	G318-2013
	2014	30	20	70	–	Xam318	G318-2014
	2013	30	20	70	–	Xam681	G681-2013
	2014	30	20	70	–	Xam681	G681-2014

(IDEAM. .

The inoculation was conducted by puncturing the stem of each plant between the second and third true leaf. The bacterial suspension was placed using a tip filled with 10 μl of the inoculum (1 × 106 UFC/mL). As mock, one plant was inoculated with 10 mM of MgCl2. The disease severity was scored at 7, 14, 21, and 30 days after inoculation, using a rating from 0 to 5 according to the symptoms scale proposed by Verdier et al. (1994), where 0 = no symptoms, 1= necrosis at the inoculation point, 2= stem exudates, 3 = one or two wilted leaves, 4 = more than three wilted leaves and 5 = plant death. Disease progress in time was estimated for each replicate by calculating the area under disease progress curve (AUDPC) (Shaner and Finney, 1977; Jeger and Viljanen-Rollinson, 2001). Once the phenotyping evaluation was finished, the inoculated material and the substrate (soil) were burned to ensure the CBB-free field condition. In order to determine resistance and susceptibility to the strain in parents and F1 genotypes, the AUDPC value was taken into account as well as the criteria previously described (Jorge et al., 2000; Restrepo et al., 2000; Trujillo et al., 2014). Even though cassava resistance to CBB is quantitative, there are several reports how to define a response as resistant or susceptible. According to Restrepo et al. (2000), if one of five biological replicate shows a scale of symptoms of 4 or 5 at 30 dpi, the genotype is considered as susceptible. On the other hand, Jorge et al. (2000) classify a genotype as susceptible when the average of the symptom using the scale, is greater than or equal to 3. Finally, the criteria described by Trujillo et al. (2014), the cassava genotypes with logAUDPC values lower than 1.59 belong to resistant genotypes while those with AUDPC values higher than 1.69 to susceptible ones.

Statistical analysis

For each location, including greenhouse, year and strain, five biological replications (clonal plants) per genotype and mock were disposed according to a randomized complete design. A Log transformation for the AUDPC values were done as described by Restrepo et al. (2000). The normality of AUDPC values were tested using Shapiro-Wilkinson test, whereas the correlation coefficients for the phenotypic data between different environments were calculated using the Pearson method. The variance components were estimated using ANOVA. All the above statistical analyses were performed using the R package.

Broad-sense heritability (H2) of the response to Xam strains was determined by calculating the genetic variance component and the residual or environmental variance component, using lme4 package (Bates et al., 2015) in R package (Ripley, 2001). The variance components of genotype, genotype by environment, genotype by year and experimental error, were measured according to Holland (2006). The value of heterobeltiosis (Jinks and Jones, 1958), was measured in order to establish the percentage of progeny that exhibited higher levels of resistance than the resistant parent. In order to evaluate performance (response to bacteria) of cassava genotypes under CBB incidence, a genotype x environment analysis was performed through a GGE-biplot analysis of multi-environment trials (MET) separately by bacteria strain and year of evaluation (rainy and dry seasons), based on the model for two principal components (Yan et al., 2000; Yan and Tinker, 2006). Data from all locations, including greenhouse, were combined to construct the GGE-biplots in order to compare the behavior of the genotypes. The GGE-biplot criteria were centered by two (centered G + GE), without scaling and singular value partitioning (SVP). GGE-biplot analysis was performed using the R package GGEBiplotGUI (Frutos-Bernal and Galindo, 2012).

QTL mapping

The QTL mapping analysis was performed with Composite Interval Mapping (CIM) using the Haley–Knott regression approach in R/qtl V1.37-11 (Broman, 2015) by employing three markers as covariate in a window size of 10 cM. The genetic positions of the 2,571 polymorphic heterozygous SNP markers are based on the genetic map distributed on 19 linkage groups with an overall size of 2,571 cM and an average distance of 1.26 cM between markers (Soto et al., 2015). A LOD threshold was calculated using 1,000 permutation tests. The LOD peak of marker-traits associations was taken into account to define the QTL position on the linkage groups in the map. The QTL interval was confined by a LOD decrease value of 1.5 compared to the peak. The phenotypic variation (R2) explained by each QTL was determined through the function calc.Rsq in R package. Environmental effects to the QTL were estimated by using the significant (LOD 2,5) additive phenotypic effects (APE) of QTL and performed using the software QTL IciMapping (Meng et al., 2015).

Identification of candidate genes

Candidate genes were located using the information of physical positions of the SNP-based genetic map. BLAST (Basic Local Alignment Search Tool) analysis was conducted to the current cassava genome sequence draft (v6.1) implemented at the JGI's Phytozome platform.

Results

Evaluation of CBB resistance in mapping population

The F1 population and parents were evaluated in two different agro-ecological locations (Arauca and La Vega) during two consecutive dry and two rainy seasons. All AUDPC values for the F1 genotypes for each location, Xam strain and season showed a continuous and normal distribution, except for LV681-D (Supplementary Material 1). As was expected, the parents exhibited AUDPC values in both extremes of the distribution curve. At the end of the trial the genotype TMS30572 was considered as resistant with AUDPC values ranging from 1.21 for LV318-R, to 1.39 for AR318-R. On the other hand, the susceptible genotype CM2177-2, exhibiting high AUDPC values of 1.72 for AR318-R to 1.92 for LV681-R. The AUDPC values observed among the genotypes were found to be highest in AR318-D (ranged from 0.97 to 1.78), followed by G681-2013 (from 1.08 to 1.89) and AR318-R (from 1.23 to 1.78) in that order. The number of genotypes evaluated, the AUDPC of the parents and the distribution of AUDPC values in the progeny for each location are shown in Table 2.

Table 2

Distribution of AUDPC values in the mapping population.

Location	Genotypes evaluated	AUDPC CM2177-2	AUDPC TMS30572	AUDPC range	Resistant genotypes	Susceptible genotypes	Resistant transgressive phenotype	Susceptible transgressive phenotype
AR318-R	103	1.72	1.39	1.23–1.78	79	24	25	4
AR318-D	100	1.76	1.26	0.97–1.78	84	16	11	1
AR681-R	104	1.83	1.35	1.28–1.94	68	36	6	2
AR681-D	100	1.8	1.32	1.16–1.94	72	28	8	3
LV318-R	93	1.85	1.21	1.16–1.82	67	26	3	0
LV318-D	106	1.73	1.24	1.04–1.82	80	26	15	7
LV681-R	93	1.92	1.39	1.19–1.97	55	38	6	2
LV681-D	106	1.88	1.36	1.21–1.97	74	32	12	6
G318-2013	117	1.85	1.28	1.09–1.88	78	39	15	1
G318-2014	109	1.87	1.28	1.10–1.86	73	36	16	0
G681-2013	117	1.87	1.25	1.08–1.89	74	43	18	1
G681-2014	112	1.87	1.25	1.10–1.90	73	39	18	2

Number of genotypes evaluated, AUDPC values for parents, range of AUDPC values obtained in each location, number of resistant and susceptible genotypes in the population for each location and number of resistant and susceptible transgressive phenotypes by environment. AR, Arauca; LV, La Vega; G, greenhouse; R rainy season; D, dry season; 318, Xam318; 681, Xam681.

In all the locations and seasons tested, the number of resistant genotypes was higher than that of the susceptible ones for both Xam strains (Table 2). Results indicated that more genotypes exhibited resistance to CBB during the dry season compared with the rainy season. About 84 and 75% of the progeny were resistant to Xam318 under dry season in Arauca and La Vega, respectively; but were reduced to 77% in Arauca and 72% in La Vega during the rainy season. In dry season, 72 and 70% of the progeny showed resistance against Xam681 in Arauca and La Vega, respectively. These values dropped to 65 and 59%, respectively, during the rainy season (Table 2).

In order to detect possible differential response of the genotypes against the pathogen, the phenotypic plasticity was determined. In at least one location and/or season, 48% of the genotypes showed a distinct phenotype for the strain Xam318, while 38% were observed for Xam681 (Supplementary Material 2). For instance, genotypes g5, g40, g53, and g116 were susceptible in AR318-R but were resistant under all the other conditions tested. On the other hand, genotypes g15, g51, and g97, which showed a resistant response to Xam681 in La Vega during rainy and dry seasons, were susceptible under the rainy season in Arauca.

The mapping population exhibited a transgressive segregation for resistance to Xam strains. Transgressive genotypes with higher resistance or susceptibility than the parents were identified in the two locations, against the two Xam strains and also under rainy and dry seasons (Table 2). However, no statistically significantly differences between the transgressive segregants and its parents were found. Locations LV318-D and LV681-D provided the highest number of transgresses, with seven and six respectively. For all conditions, the total number of resistant transgresses was higher than susceptible ones. The highest number of resistant transgresses was identified at AR318-R with 25 transgresses and in the greenhouse conditions against Xam681 in both G681-2013 and G681-2014, with 18 and 18 transgresses each (Table 2). The genotype g79 exhibited a resistant transgressive phenotype in six of the 12 conditions evaluated (Arauca, for both Xam strains, and also during rainy and dry seasons). The greenhouse evaluation using Xam681 showed lower AUDPC values in relation to the resistant parent. However, the genotype g29 was identified as a susceptible transgressive phenotype for six of the 12 conditions: AR681-R, AR681-D, LV681-D, G318-2013, G681-2013 and G681-2014 (Supplementary Material 2).

Pairwise Pearson correlation between AUDPC values measured at two locations and greenhouse, two seasons and for the two Xam strains were highly significant (P < 0.05) with correlation coefficients ranging between 0.62 for pairwise comparison of AR318-R and AR318-D to 0.99 for G681-2013 and G681-2014 (Table 3). Nevertheless, very low values of correlation were found between different environments.

Table 3

Pairwise Pearson correlation coefficients between AUDPC values for all environments and for the two Xam strains employed.

Location	AR318-R	AR681-R	LV318-R	LV681-R	G318-2013	G681-2013
AR318-D	0.62
AR681-D	0.06	0.84
LV318-D	0.03	0.08	0.81
LV681-D	0.15	0.35	0.05	0.89
G318-2014	0.07	0.15	0.04	0.05	0.91
G681-2014	−003	0.34	0.04	0.22	0.34	0.99

AR, Arauca; LV, La Vega; G, greenhouse; R, rainy season; D, dry season; P-value = 0.05

Broad sense heritability of the resistance to Xam318 was 20%, with values of 0.0022 and 0.0110 for genetic and environmental variance, respectively. While the broad sense heritability of the resistance to Xam681 was 53%, with values of 0.0075 and 0.0141 for genetic and environmental variance, respectively. The highest and most stable values of heterobeltiosis were observed for the tests conducted under greenhouse conditions. For strain Xam681 the values were 21.11 and 20.39% in 2013 and 2014, respectively, whereas for Xam318 it was 16.62% for both years of evaluation (Table 4).

Table 4

Better-parent heterosis.

	AR-R (%)	AR-D (%)	LV-R (%)	LV-D (%)	G-2013 (%)	G-2014 (%)
Xam681	14.36	14.60	13.46	14.15	21.11	20.39
Xam318	7.43	14.24	14.60	16.59	16.62	16.62

Percentage of heterosis over the better parent values for the quantitative genetic trait of resistance to the two Xam strains evaluated in cassava under the environments AR, Arauca; LV, La Vega; G, greenhouse; R, rainy season; D, dry season.

The analysis of variance showed significant differences (p < 0.001) for genotype (g), environment (location, season), Xam strain, and genotype x environment among the F1 genotypes tested (Supplementary Material 3). The GGE-Biplot analysis revealed that the two principal components explained 77.83% of the total variance caused by Gplant + Gpathogen + E for cassava resistance to Xam318 during rainy season, 83.96% for Xam318 during dry season, 81.12% for Xam681 in rainy season and 83.97% for Xam681 in dry season (Figure 1). Also, the AUDPC values of the genotypes were able to discriminate between environments. The behavior of genotypes differed between environments and these fall into two sectors of the graphic, except for Xam318 during the rainy season, which fall in three sectors (Figure 1). In both Arauca and La Vega, except for Xam318 rainy season, the genotypes behavior seems to be quite similar (Figures 1B–D).

Figure 1

Which Won Where/What graphic of GGE-Biplot analysis. (A) GGE-Biplot analysis for phenotype evaluation against Xam318 under rainy season. (B) GGE-Biplot analysis for phenotype evaluation against Xam318 under dry season. (C) GGE-Biplot analysis for phenotype evaluation against Xam681 during the rainy season. (D) GGE-Biplot analysis for phenotype evaluation against Xam681 during the dry season. e1, Arauca; e2, La Vega; e3, Greenhouse.

The environments with short vectors indicate that the genotypes under this environment behave similarly. This is the case for the evaluation in Arauca (e1) against Xam318 during rainy and dry season (Figures 1A,B). In the same category falls the evaluation conducted for Xam681 in Arauca during rainy season (Figure 1D). On the other side, large vectors such as those exhibited by all the evaluations under the greenhouse condition, indicates higher level of dissimilarity between the phenotypic responses of the genotypes (Figure 1).Using the criteria previously described (Jorge et al., 2000; Restrepo et al., 2000; Trujillo et al., 2014), it was possible to distinguish genotypes as the “extremes” for resistant and susceptible for almost all environments. Genotypes g33 and g41 were selected as extreme resistant and susceptible, respectively in response to Xam318 in rainy season. During the dry season, the extreme resistant genotype against the same strain was g131, while g124 was identified as the extreme susceptible genotype. On the other hand, the evaluation conducted under the rainy season employing Xam681 allowed to distinguish g79 and g29 as extreme genotypes for resistance and susceptibility, respectively. The genotype g93 was identified as the extreme resistant and g29 as the extreme susceptible in the dry season. The g29 showed to be the most susceptible for the two Xam strains and it also had the highest susceptible transgressive phenotype.

Identification of QTL for resistance to Xam in cassava

In total 5 QTL distributed in 4 of the 19 linkage groups (LG) were detected (Figure 2). The QGH681-2.2 located on the LG 2.2, the QLV681RD-6 on the LG 6, the QAR681D-14 on LG14 and the QLV318RD-19 and QGH318-19 located on the LG 19. The QTL detected in the evaluation at La Vega during the rainy season (QLV681RD-6) explaining 22.1% of the resistance to Xam681, showed the highest LOD value of 5.0 (Table 5). The phenotypic variance of resistance to Xam explained by all identified QTL ranged from 15.8 to 22.1%. The evaluation against Xam318 detected two QTL, explaining 17.3–18.8% of the phenotypic variance, whereas three QTL were identified for resistance to Xam681 accounting for 15.8–22.1% of the phenotypic variance. From the 5 QTL, one was detected at Arauca location, whereas two QTL were identified each at Vega and greenhouse conditions.

Figure 2

Schematic view of the QTL for resistance to CBB located within cassava linkage groups. In the represented cassava linkage groups are shown the detected QTL, in red letters are shown the QTL peak markers, in black, the flanking markers. The genetic distances in cM are shown with the scale on the left. Diagram plotted using MapChart software (Voorrips, 2002).

Table 5

Summary of QTL for resistance to CBB.

QTL name	LOD score	R2	Peak Marker	Pos. cM	Marker Interval	Interval length cM	Interval length Kb	Genes in interval
QAR681D-14	4.3	18.1%	MB_9956	90.6	MB_62539/MB_10857	4.1	7.1	2
QLV318RD-19	3.8	17.3%	MB_38006	*	MB_37785/MB_38030	*	80.8	8
QLV681RD-6	5.0	22.1%	MB_50599	31.0	MB_36352/MB_63835	2.2	6.5	3
QGH318-19	4.9	18.8%	MB_49863	*	MB_49851/MB_49878	*	58.8	8
QGH681-2.2	4.2	15.8%	MB_23160	93.5	MB_23143/MB_26945	2.3	110.1	8

QTL name, location, Xam strain and linkage group), LOD threshold, R.

Unknown genetic position.

On average, the interval length of the QTL was 2.8 cM, with the lowest interval length 2.2 cM for QLV681RD-6 and the largest of 4.1 cM for QAR681D-14 (Table 5). The QTL QGH681-2.2, QLV681RD-6, QLV318RD-19, and QGH318-19 were stable between the two seasons for the same location. From these, two correspond to QTL detected under greenhouse conditions and two under La Vega conditions. Two of the stable QTL confer resistance to Xam681 and the other two to Xam318 (Table 5). There were no QTL detected for all environments and/or seasons with LOD score above the threshold.

Interaction of QTL with environment

In order to estimate the environment effects on the QTL, the Q × E interaction analysis was conducted. In this case, the APE obtained through the QTL IciMapping was emloyed as a measure of the environment effect on the QTL identification. The QTL detected under greenhouse conditions were not taken into account. Three QTL exhibit significant differences of APE between the evaluated environments, these environment specific QTL effects were interpreted as Q × E interactions (Figure 3).

Figure 3

QTL x environment interactions based on additive phenotypic effects (APE). QTL exhibit significant Q x E interactions when considering significant (LOD 2,5) additive phenotypic effects (APE). Positive APE was correlated with susceptibility to Xam while negative APE with resistance.

The QLV318RD-19 exhibits positive APE, which was correlated with susceptibility effect or high AUDPC values from the phenotypic evaluation. However, the QLV681RD-6 and QAR681D-14 showed negative APEs and correlate well with resistance or low AUDPC values. The QAR681D-14 is a conditionally neutral QTL, because it was detected only in a specific environment. The QTL QLV681RD-6 and QLV318RD-19 are stable, showing different effect levels in dry and rainy seasons. However, the QTL QLV681RD-6, detected during rainy and dry seasons, showed a significant APE value only under dry conditions (Figure 3).

Sequence analysis of chromosome regions containing the identified QTL

The DNA sequences of genome intervals harboring each QTL were analyzed in silico for the presence of genes by Blasting in the public databases. The found hits were considered as CBB CDRGs. In total, 29 CDRGs were found and they are distributed in all the regions containing QTL (Supplementary Material 4). The QTL physical intervals comprised 263.3 Kb corresponding to a genetic distance of 8.6 cM with a mean value of 30.6 Kb per 1 cM. On average, one gene each 9.07 Kb was found. However, this ratio varies between the QTL from 0.07 in QGH681-2.2 region to 0.4 genes per Kb in QLV681RD-6. QTL intervals QGH681-2.2, QLV318RD-19, and QGH318-19 which showed the highest interval-lengths (101.1, 80.0, and 58.8 Kb, respectively) contained the highest number of genes (8 genes in each interval). However, the gene density (GD) varied from one gene every 7.3 Kb for QGH318-19 region to 13.7 Kb for QGH681-2.2.

Among the 29 CDRGs co-localizing with the QTL, five genes (17.2%) have not been previously annotated and/or described in the current cassava genome, based on PFAM (Finn et al., 2014), PANTHER (Thomas et al., 2003) and EuKaryotic Orthologous Groups (KOG) (Koonin et al., 2004). The identified candidate genes contain diverse protein domains, without particular domains being most represented. Interestingly, two immunity related genes (IRGs) co-localized with QLV318RD-19, but 1.4 Mb apart from each other: one of them coding for a WRKY DNA -binding domain bearing protein and the second coding for a NB-ARC-LRR class of proteins. A gene coding for a subtilisin-like protease was found within the interval of the stable QTL QLV681RD-6. Also, a gene containing kinase domain was identified in the QGH681-2.2 interval (Supplementary Material 4).

Discussion

In this study, the defense response of 117 genotypes of a F1 mapping population infected with two Xam strains was evaluated in different environments: two locations under rainy and dry seasons and under greenhouse controlled conditions, which allowed for the identification of five strain-specific QTL for resistance to CBB. In silico analyses of the QTL intervals revealed 29 CDRGs that might operate in the defense response. The strains used in this study are part of a set of Xam strains evaluated for virulence in nine cassava accessions. The Xam318 and Xam681 showed a high virulent behavior, causing disease in six and eight of nine accessions tested, respectively (Trujillo et al., 2014). Moreover, these strains have been reported as part of prevalent haplotypes described in Colombian Xam populations and to be the product of migratory processes between regions in Colombia (Trujillo et al., 2014). The Xam populations in these regions remained unexplored since more than one decade (Restrepo, 1999). Although cassava is cultivated in both regions, the grown areas are relatively limited. It remains interesting to examine if the two Xam strains used for the resistance test in this study do exist in the Xam populations in La Vega, Cundinamarca and Arauca where the studies were performed. This will allow directing the breeding programs toward the resistance to these Xam strains for these particular Colombian regions.

The broad sense heritability of cassava plants estimated in this study was 20 and 53% for resistance to Xam318 and Xam681, respectively. The values are in the range of the previously reported studies of 10–69% (Hahn et al., 1990; Jorge et al., 2000; Fregene et al., 2001; Ly et al., 2013). These values of heritability highlight the important effect of the environmental conditions on plant response to CBB. In particular, the humidity showed a strong influence on the phenotypic response of some of the F1 individuals. Thus, the number of susceptible individuals was higher during the rainy season compared to the dry season in both environments (68 vs. 64 for LV and 60 vs. 44 for Arauca). Several studies have shown a positive (favorable) effect of humidity not only on the speed of symptoms, but also on the growing of Xam (Banito et al., 2001; Wydra and Verdier, 2002; Restrepo et al., 2004). Even though the environments with higher humidity seem to favor the disease, temperature also can generate different effects under a particular genotype. This can lead to a change in the plant phenotype, which is defined as phenotypic plasticity (Nicotra et al., 2010). In all the conditions evaluated genotypes that exhibited resistant behavior under some environmental conditions but susceptible in others were found. For example, three genotypes g29, g30, and g116, showed resistance under the dry season but susceptibility under the rainy one, whereas g51 and g97 which were resistant to Xam681 in La Vega during rainy and dry seasons, but susceptible under the rainy season in Arauca. The phenotypic plasticity has been widely described in model and non-model crops for several traits, including resistance to plant pathogens (Agrawal, 1999; Dicke and Hilker, 2003). The individuals showing phenotypic plasticity identified in this study can be used in local breeding programs as superior genotypes adapted to specific environmental conditions as an approach exploiting adaptive plasticity (Nicotra et al., 2010).

Despite that no statistically significantly differences between the transgressive segregants and its parents were found, it is important to consider a clear differential biological behavior of CBB resistance and symptoms for some of these genotypes. In this way and taking into account other factors, these genotypes represent important sources of resistance for future use in cassava breeding programs. In addition, these genotypes can be also considered for future research aiming to better understand the mechanisms of bacterial resistance in cassava. Transgressive segregation of resistance has been previously described in cassava against Xam strains (Jorge et al., 2000) as well as for other important traits (Akinwale et al., 2010; Whankaew et al., 2011; Njenga et al., 2014; Thanyasiriwat et al., 2014). This type of heterosis has been explained by the presence of blocks of dominant genes from both parents (Bingham, 1998; Jorge et al., 2000), variations in the chromosome number, chromosome rearrangements (Rieseberg et al., 1999), and even DNA methylation, epigenetics and silencing by small RNAs (Shivaprasad et al., 2012). The highest resistance was scored in the evaluations conducted under controlled conditions. A particular case was the genotype g79, which was categorized as a resistant transgressive phenotype for almost all the conditions evaluated. The genotypes g3, g103, g118, g126, were detected showing extreme resistance in both rainy and dry season against Xam318 and Xam681 (Supplementary Material 2). These genotypes could be important sources of resistance for future use in cassava breeding programs.

Five strain specific QTL conferring resistance against Xam in cassava were identified in this study. Three of the associated QTL were found to be effective against Xam318 strain, while two were detected for Xam681. Jorge et al. (2000), have reported 12 QTL for resistance to five Colombian Xam strains, while Wydra et al. (2004) nine QTL to four Xam strains. None of these QTL were identified in the present study, and could be explained by: the lack of common markers between the used genetic maps; by the very different environments where the evaluations were conducted; by the different genetic backgrounds of the mapping populations. The strains used in this study are different from those used in above mentioned reports, and all the QTL we identified were strain-specific. Detection of strain-specific QTL has already been reported for quantitative resistance in several crops such as rice (Li et al., 1999), tomato (Wang et al., 2000; Carmeille et al., 2006), melon (Perchepied et al., 2005), and apple (Calenge et al., 2004). Taken together, these results suggest a strong strain x cultivar x environment interaction. In order to prove this hypothesis, it will be important to consider phenotypic evaluations in multi-environments with the same Xam strain as well as to expand the repertoire of strains belonging to the same and different haplotypes.

The use of the APE as a parameter to establish Q x E interaction have been reported in some crops such as wheat (Hao et al., 2011) and in the plant model Arabidopsis (El-Soda et al., 2014). Here, the Q × E interactions were established for all QTL except for those detected under greenhouse conditions. Three QTL showed significant APE values. Remarkably, the QTL detected under dry conditions exhibit a negative APE. This shows once again the important role of the environment and especially of the humidity in favoring CBB disease. The detection of QTL with a significant APE in one environment but not in another, unstable QTL or conditionally neutral QTL, are evidence of the pivotal role of the environment conditions in the instability of the detection of QTL between environments. Conditionally neutral resistance QTL have been reported in some crops like rice (Li et al., 2007), wheat (Ramburan et al., 2004), and apple (Calenge and Durel, 2006). In cassava, in spite of that environmentally unstable QTL for some traits have been reported (Jorge et al., 2001), Q × E interaction for CBB has not been described so far. Although, stable detected QTL are most useful for breeding programs, the mechanisms behind a conditional neutral QTL are useful in plant improvement. These resistance loci indicate the great influence of external conditions, such as environmental factors on the phenotype, thus, they can be exploited in local programs that present particular environmental conditions and adapted pathogens.

All the QTL identified cover in total a physical part of the cassava genome corresponding to 263.3 Kb. In these regions 29 CDRGs are found. With the advent of new high throughput genotyping technologies the number of markers increased exponentially and allowed the construction of genetic high resolutions maps in cassava (Soto et al., 2015). The efforts made recently to sequence the cassava genome allow a better characterization of the QTL intervals, including the number and nature of the genes present in these regions. In spite of associations between candidates genes with QTL have been reported (Faris et al., 1999; Ramalingam et al., 2003), the isolation of genes from QTL regions are scarce. By using a genetic map with higher resolution we were able to identify 29 CDRGs that genetically co-localized with QTL in short interval lengths (2.8 cM in average).

From the 24 annotated genes, two are related to plant immunity and showed to be co-localized with one QTL (QLV318RD-19). Similar findings have been reported for other plant species (Ramalingam et al., 2003; St. Clair, 2010; Lopez, 2011). The resistance to CBB has been classified as a quantitative trait and in cassava Avr-R interactions have not been demonstrated so far. Nevertheless, the presence of genes coding for typical R proteins within the QTL intervals supports the finding of an overlapping between qualitative and quantitative resistances.

In the whole repertoire of genes present in the QTL intervals, a kinase co-localized with QGH681-2.2. Despite the kinase family being one of the most widely distributed proteins in plant genomes, finding a co-localization of a member of this family with resistance seem not to be by chance. Several proteins of this family have been shown to be involved in plant resistance in model plants (Huard-Chauveau et al., 2013) as well as in some of the most important crops such as barley (Druka et al., 2008), wheat (Fu et al., 2009), and maize (Zuo et al., 2015). The kinases can be an important element in quantitative disease resistance, either as a receptor or as part of the signaling pathway. Kinase proteins have been reported to be involved in defense response pathways in Arabidopsis and in other crops (Roux et al., 2014). In tomato for example, the Pseudomonas resistance gene (Pto) encodes a Ser/Thr kinase (Oh and Martin, 2011). The large group of kinases has been exploited as “decoys” for the study of plant immunity (van der Hoorn and Kamoun, 2008).

Even though some of the CDRGs present in the QTL regions do not belong to IRGs, they still need to be considered as potential candidate genes for CBB resistance. Conversely, it has been demonstrated that the classical immunity-related proteins could be only a small fraction of the total genes involved in quantitative resistance (Corwin et al., 2016). Cloning of some genes for plant quantitative resistance has given evidence that the corresponding proteins do not belong to any specific group of immunity related proteins or they lack the conserved domains of R proteins. These genes are involved in different functions and processes (Poland et al., 2009; Bryant et al., 2014; Roux et al., 2014). Thus, the gene coding for a subtilisin-like protease that co-localized with the stable and major QTL QLV681RD-6, could be considered as an interesting and promising CDRG for several reasons: first, due to fact that it co-localizes with a QTL that explains the highest percentage of phenotypic variance; second, the coding protein subtilisin-like protease has been described to be involved in resistance to pathogens in other crops such as tomato (Tornero et al., 1996; Jordá et al., 1999), grapevine (Figueiredo et al., 2008; Monteiro et al., 2013) and in Arabidopsis (Ramirez et al., 2013; Figueiredo et al., 2014); third, this candidate gene co-localizes with a stable QTL, suggesting that the effect in conferring defense response to Xam of this gene is not strongly affected by the environment.

The repertoire of CDRGs co-localizing with the QTL reported here represents a first step in the dissection of the molecular mechanisms that govern CBB resistance in cassava and a new source of genes to be validated through different approaches in the future. With the advent of gene editing methodologies it will be possible to validate the function of these particular genes in CBB resistance (Sander and Joung, 2014). Moreover, this gene repertoire and the associated SNP markers can be used in plant breeding to develop CBB resistant cassava that are adapted to the evaluated regions and other locations with similar environmental conditions.

Author contributions

JS carried out the CBB phenotyping and QTL analysis. RM conducted part of the CBB phenotyping. JS, AB, and CL wrote the manuscript and contributed to the interpretation of the results. BM and JL performed the statistical analysis. JS and CL designed the experiments. CL coordinated the project. FG performed the CDRGs analysis. All authors read and approved the final manuscript.

Conflict of interest statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.