Email updates

Keep up to date with the latest news and content from BMC Plant Biology and BioMed Central.

Open Access Research article

Integration of promoters, inverted repeat sequences and proteomic data into a model for high silencing efficiency of coeliac disease related gliadins in bread wheat

Fernando Pistón*, Javier Gil-Humanes and Francisco Barro

Author Affiliations

Instituto de Agricultura Sostenible, Consejo Superior de Investigaciones Científicas (IAS-CSIC), Córdoba E-14080, Spain

For all author emails, please log on.

BMC Plant Biology 2013, 13:136  doi:10.1186/1471-2229-13-136

The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1471-2229/13/136


Received:12 March 2013
Accepted:4 September 2013
Published:17 September 2013

© 2013 Pistón et al.; licensee BioMed Central Ltd.

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

Background

Wheat gluten has unique nutritional and technological characteristics, but is also a major trigger of allergies and intolerances. One of the most severe diseases caused by gluten is coeliac disease. The peptides produced in the digestive tract by the incomplete digestion of gluten proteins trigger the disease. The majority of the epitopes responsible reside in the gliadin fraction of gluten. The location of the multiple gliadin genes in blocks has to date complicated their elimination by classical breeding techniques or by the use of biotechnological tools.

As an approach to silence multiple gliadin genes we have produced 38 transgenic lines of bread wheat containing combinations of two endosperm-specific promoters and three different inverted repeat sequences to silence three fractions of gliadins by RNA interference.

Results

The effects of the RNA interference constructs on the content of the gluten proteins, total protein and starch, thousand seed weights and SDSS quality tests of flour were analyzed in these transgenic lines in two consecutive years. The characteristics of the inverted repeat sequences were the main factor that determined the efficiency of silencing. The promoter used had less influence on silencing, although a synergy in silencing efficiency was observed when the two promoters were used simultaneously. Genotype and the environment also influenced silencing efficiency.

Conclusions

We conclude that to obtain wheat lines with an optimum reduction of toxic gluten epitopes one needs to take into account the factors of inverted repeat sequences design, promoter choice and also the wheat background used.

Keywords:
Silencing; Endosperm-specific; RNAi; Coeliac; Wheat; Gluten

Background

Wheat is one of the most important food grains in the world, being processed into bread and many other products. Wheat products not only make substantial contributions to the dietary intake of energy and protein, but also have impacts on human health, both beneficial (providing dietary fibre, minerals, vitamins, phytochemicals) and negative (allergies and intolerances, which appear to be increasing in importance). Of particular interest is the gluten fraction of wheat grain, as this not only imparts unique technological characteristics, but also plays a major role in intolerance and allergy. The gluten proteins represent 80% of the total grain protein in bread wheat (Shewry and Halford, 2002), and comprise two major groups: the glutenins and the gliadins [1]. The glutenins include the high molecular weight (HMW) and the low molecular weight (LMW) fractions, whereas the gliadins can be divided into three structural types: α-, ω-, and γ-gliadins [2].

In recent years, there have been increases in the incidence of wheat allergies and intolerances. Food allergy to wheat affects 0.2-0.5% of the population [3] and is difficult to manage because very many food products contain wheat. Wheat allergies affects adults and children, can cause anaphylaxis and are triggered by the consumption of gluten. Coeliac disease (CD), which is an intolerance to glutens from wheat, rye and barley and occurs in both children and adults throughout the western world at an average frequency of about 1%, with some groups reporting rates five-fold greater in infants compared with adults [4]. Gluten sensitivity is a new pathology of intolerance to gluten [5], which excludes CD and wheat allergy. No accurate estimates on the prevalence of gluten sensitivity are available, but preliminary data for the USA (6% of the population) suggest that it is more frequent in the general population than CD. CD and wheat allergies require a strict gluten-free diet and sufferers from gluten insensitivity aim to reduce gluten intake as far as possible. Therefore, the development of novel gluten-free wheat varieties is a major objective. CD is the most studied of the gluten-related pathologies. It is an autoimmune disorder with genetic and immunological components as a consequence of the ingestion of gluten proteins from wheat and related cereals. The peptides produced in the digestive tract by the incomplete digestion of gluten cause inflammation of the small intestine and villous atrophy. The autoimmune response is a consequence of the deamidation of glutamine residues present in the peptides, by the tissue transglutaminase 2 (tTG2) in the gut mucosa. The deaminated peptides are able to bind to class II human histocompatibility leukocyte antigen (HLA) molecules DQ2 and DQ8, which stimulate T cells, resulting in an inflammatory response in the small intestine that leads to flattening of the mucosa [6]. Isolation and characterization of intestinal T cells from CD patients have revealed several distinct but similar DQ2 and DQ8 epitopes. Although a number of epitopes are derived from glutenins [7], the majority of the epitopes reside in the gliadin fraction [8,9].

Wheat gliadin genes occur in tightly linked clusters, termed blocks, located at complex loci on group 1 and 6 chromosomes [10]. The estimated copy numbers in hexaploid wheat of genes encoding α-gliadins ranges from 25 to 150 copies [11], from 15 to 18 copies for ω-gliadins, and from 17 to 39 copies for γ-gliadins [12]. This high level of complexity [11,12] and the fact that gliadin genes are inherited in blocks make conventional breeding approaches to obtain wheat varieties with reduced content of T-cell stimulatory epitopes very unlikely.

Post-transcriptional gene silencing by RNA interference (RNAi) is based on sequence-dependent RNA degradation that is triggered by the formation of double-stranded RNA (dsRNA), homologous in sequence to the targeted gene [13]. We used this approach to down-regulate the expression of coeliac disease-related wheat gliadin T-cell epitopes [14,15], showing that RNAi technology can be used to obtain wheat varieties free of toxic epitopes and potentially suitable for CD patients. The use of specific inverted repeat sequences and optimal endosperm-specific promoters, which drive the expression of hairpin constructs, are critical factors in achieving effective down-regulation of CD-related gluten proteins and to minimize off-target effects.

In the present study, 38 transgenic wheat lines with different gliadin fractions down-regulated by RNAi were analyzed over two consecutive years The effect of the promoters and RNAi fragments used for silencing, the genotype, as well as the effect of environment on the content of the gluten protein fractions, total protein and starch, thousand seeds weight and SDSS quality test of flour were analyzed.

Results

In this work we report the assay of 38 transgenic lines with different gliadin fractions down-regulated by RNAi using two different endosperm specific promoters (D-hordein and gammma-gliadins promoter) and two silencing fragments (to target gamma-gliadins and all gliadins fractions) in a randomized complete block design (RCBD). We used the four combinations of promoter and silencing fragment. The following constructs were used: the vector pghpg8.1 has the D-hordein promoter and the silencing fragment against gamma-gliadins; the vector pGghpg8.1 has the gamma-gliadin promoter and the silencing fragment against gamma-gliadins; pDhp-ω/α vector has the D-hordein promoter and the silencing fragment to target all gliadins fractions; the pGhp-ω/α vector drives the silencing fragment to down-regulate all gliadin fractions with the gamma-gliadin promoter. The assay was repeated during two years to evaluate the environment interaction and gluten protein composition, the thousand seed weights, total protein and starch, and SDSS test were evaluated in transgenic and control lines.

Analysis of protein fractions from gluten

The gliadin and glutenin fractions from the 38 transgenic lines and their respective wild type controls were analyzed by RP-HPLC, where proteins are eluted according to different surface hydrophobicity. The elution order is ω-, α-, and γ-type for the gliadin fraction and HMW and LMW subunits for the glutenin fraction [16,17]. Figure 1 shows the pattern obtained for BW208 and BW2003 controls and four transgenic lines with gliadins down-regulated (Figure 1a and b). The control lines BW208 and BW2003 had different patterns for the gliadin and LMW-GS fractions but the same HMW-GS composition (Figure 1c and d). Chromatograms from transgenic lines were compared with their respective control lines, showing clear differences. Lines A1158 and 22A, transformed with the pghpg8.1 and pGghpg8.1 vectors targeting the γ-gliadins, showed a strong decrease of the γ-gliadin region, whereas the ω- and α-gliadin regions increased the areas of their peaks in the transgenic lines of both genotypes (Figure 1a and b). In the glutenin fractions, these transgenic lines showed an increase of the LMW-GS and HMW-GS (Figure 1c and d). Lines E33, transformed with the pDhp-ω/α and pDhpg8.1 vectors, and D874, transformed with the pDhp-ω/α and pGhp-ω/α vectors, showed a reduction of the α-, ω-, and γ-gliadins region. Both lines E33 and D874 showed an increase in the HMW-GS but line D874 also showed reduction in the LMW-GS fraction (Figure 1d).

thumbnailFigure 1. RP-HPLC chromatograms of gliadin and glutenin extract representing wild-type control and transgenic wheat lines. (a) Gliadin extracts from BW208 wild-type and BW208 transgenic lines. (b) Gliadin extracts from BW2003 wild type and BW2003 transgenic lines. (c) Glutenin extracts from BW208 wild type and BW208 transgenic lines. (d) Glutenin extracts from BW2003 wild type and BW2003 transgenic lines. ω, ω-gliadins; α, α-gliadins; γ, γ- gliadins; HMW, high molecular weight glutenin subunit; LMW, low molecular weight glutenin subunit. mAU (210 nm), milliunits of absorbance at 210 nm; min, retention time in minutes.

We carried out the quantitative analysis of gliadin and glutenin content of all 38 transgenic and wild type control lines by the integration of chromatograms. The efficiency of the constructs on gliadin silencing and their effects on gluten protein was analyzed by a linear mixed model where total protein content (measured by Kjeldahl method) was used as fix covariate, with iTarget and Genotype as fixed-effect factors and Year and Block as random effects (Model 1; Table 1). The p-value and explained variance from each effect in the model is reported in Table 1. Total protein, as covariate, had a significant effect on gliadin and glutenin fractions except on the content of LMW-GS and total glutenin content, although its contribution to the overall variance was the lowest of that of the fixed-effects. iTarget was highly significant on the content of gliadin and glutenin fractions, and it was the effect with the higher explained variance, except for the ω-gliadin content. Genotype had a significant effect on the content of ω-, α- and γ-gliadins, and on HMW-GS and LMW-GS, but not on the total content of gliadins, glutenins and Gli/Glu ratio (Table 1). The interaction between iTarget and Genotype (iTarget x Genotype) was also significant for all the gliadin and glutenin fractions, with the exception of ω-gliadins. Together, Genotype and iTarget x Genotype, had lower explained variance than the iTarget effect. The random-effect interaction iTarget x Year was significant for the ω-gliadin content. Year also had a significant effect on the ω-gliadin and LMW-GS contents and on the Gli/Glu ratio, and a borderline significance on total gliadin and HMW-GS contents (Table 1). Except for the contribution of Year to the total variance of the ω-gliadin content (69.77%), the variation of all the variables was mainly due to the fixed effects.

Table 1. Significance and percent of explained variance of the fixed and random effects of mixed model for each variable studied

Levels of factor iTarget were compared by a post hoc multiple-comparison analysis (Table 2a). The iTarget ‘g’ (constructs pghpg8.1 and pGghpg8.1) led to a decrease of the γ-gliadin content and to an increase of the HMW-GS and LMW-GS contents relative to control (iTarget C). However, total prolamin content and the Gli/Glu ratio were not significantly affected (Table 2a). The iTarget ‘o’ (constructs pDhp-ω/α and pGhp-ω/α) decreased the content of ω-, α- and γ-gliadins, total gliadins, total prolamins, and Gli/Glu ratio; and increased the content of HMW-GS and LMW-GS. The combination of constructs ‘g’ and ‘o’ (‘iTarget go’) led to a similar effect on the gliadin content than the construct ‘o’ alone, but reduced the levels of all three gliadin fractions and in particular the γ-gliadins fraction. The combination of constructs lead to a lower increment in HMW-GS content than construct ‘o’ alone, but a comparable decrease in LMW-GS content relative to control. Overall, the total glutenin content was compensated (higher HMW-GS and lower LMW-GS) and did not show any difference in comparison to control lines.

Table 2. Prolamins and non-gluten contents, gliadins and glutenins contents, total protein and starch contents, SDS sedimentation test and 1000 grain weight of transgenics grouped by silencing fragment and promoter used

To compare the effect of the two promoters used for the silencing of gliadins, the values from transgenic lines, relative to their wild type control lines, were transformed inside each Block of the experimental design (e.g. value from transgenic line of genotype BW208 of Block 1 from year 2010/control line value of genotype BW208 of Block 1 from year 2010). The fixed-effects were total protein as covariate, iTarget and Promoter; and the random-effects were Block and Year, and their interactions (Model 2; Table 1). The Promoter had a significant effect on the γ-gliadin, LMW-GS, total content of glutenins and prolamins. Furthermore, the Promoter effect explained higher variance (in comparison with the mains factors) in the γ-gliadin content, and explained also high variance in the LMW-GS and total glutenin contents.

Promoter factor was compared by a post hoc multiple-comparison of means (Table 2b). The Promoter factor has three levels, the constructs driven by the D-hordein promoter (‘D’), the γ-gliadin promoter (‘G’) and the combination of constructs with ‘G’ and ‘D’ promoters (‘GD’). As shown above, the Promoter factor had effect on the γ-gliadin, LMW-GS, total glutenin and prolamin contents. The strongest decrease in γ-gliadins was seen with the combination of promoters (‘GD’), and the least effect with the D-hordein promoter. The content of LMW-GS was higher with the D-hordein promoter and lower with the combination of both promoters. In addition, when both promoters were combined, they provided lower total glutenin content than using the γ-gliadin and D-hordein promoters independently. The prolamin content was higher with the D-hordein promoter and lower when both promoters were combined (‘GD’).

The Line effect was also evaluated by a mixed model described above (Model 3). The factor Line had a significant effect on the content of all the gliadins and glutenins. To check the differences between transgenic and control lines (each transgenic line against its control genotype), a multiple-comparison of means was made (Table 3). The comparison of means between lines showed an expected behavior according to the iTarget silencing construct used in each case (‘g’, ‘o’ or ‘go’) (see also Table 2b) for the content of ω-, α-, and γ-gliadins, total gliadins, total prolamins and the Gli/Glu ratio. Furthermore, the lines had a higher variability for the content of glutenin (HMW-GS, LMW-GS and total glutenin contents), and these did not follow the same pattern described above for each silencing construct. Thus, we found lines with the silencing constructs ‘g’, ‘o’ and ‘go’ with HMW-GS content equal to their control, in contrast to what was seen when the lines with the same iTarget were analyzed together. Similar observations were made for the content of LMW-GS and total glutenins.

Table 3. Prolamins and non-gluten contents, gliadins and glutenins contents, total protein and starch contents, SDS sedimentation test and 1000 seed weight of transgenic and wild-type lines

Thousand seed weights, total protein and starch, and SDSS assay analysis

Results from SDSS assay, thousand seed weights, total protein and total starch contents, were also analyzed using the mixed model described above (Models 1, 2 and 3; Table 1). The total protein content was analyzed as dependent variable using the same models without the covariate (Table 1). In Model 1 the iTarget had a significant effect in the four dependent variables and it was the effect with the higher explained variance. The covariate effect showed significant differences for total starch and thousand seed weights. Genotype showed significant effect on thousand seed weights and total starch. The interaction iTarget x Genotype led to significant differences in SDSS values. As for the variables analyzed previously (above), the main effect was iTarget, with the higher explained variance for all variables. However, for total starch, the Genotype and total protein together explained almost the same variance as iTarget. The more remarkable random-effects were the Genotype x Year for the SDSS and Year for thousand seed weights. A multiple comparison of adjusted means of the iTarget levels was carried out on these variables (Table 2a). Total protein content of transgenic plants transformed with the iTarget ‘go’ was significantly higher. The transgenic plants with iTarget fragments ‘o’ and ‘go’ had lower total starch contents and thousand seed weights than iTarget ‘g’ and control lines. In the SDSS analyses, the transgenic lines with iTarget ‘g’ showed an increase of SDSS values relative to controls whereas the iTarget ‘go’ showed a decrease and the iTarget ‘o’ a borderline decrease (p-value = 0.0557).

The model used to evaluate the effect of the Promoter with the data transformed by the control (Model 2) showed a significant effect of Promoter on SDSS and thousand seed weights but not on Total Starch and Total Protein contents (Table 1). No random factor showed an important effect on the variables. A multiple comparison of adjusted means of the Promoter levels was carried out on these variables (Table 2b). For the SDSS variable, the combination of both promoters was different to the other two, which showed no differences (Promoter ‘D’ and ‘G’). The promoter ‘G’ was different to the promoter ‘D’ and ‘GD’ for thousand seed weight.

Multivariate analysis of variance (MANOVA) and NMDS

To analyze the data as a whole and to determine which factors have more weight in this dataset, we conducted a nonparametric multivariate analysis of variance (MANOVA). Differences between the two Promoters were tested in the analysis with the data normalized. MANOVA results (Table 4) showed that iTarget, Year, the covariate Total protein, and the interactions Block x Year, and Year x iTarget had a significant effect on the whole data (Table 4; Model M1). The factor iTarget was the main effect and explained 42.4% of the data variability. The interaction Block x Year (5.9%), Year (5.7%), Year x iTarget (3.7%), and the covariate Total protein (1.6%) had minor effect on the data variability.

Table 4. Permutational multivariate analysis of variance comparing the whole set of plant traits between total protein, year, genotype, promoter and iTarget of silenced and control plants

The MANOVA with normalized data showed, in those common factors and covariates, similar results to the previous Promoter model (Model 2). The Promoter factor, introduced in this model, turned out to have a significant effect on the data set (Table 4; Model M2). In terms of explained variance, the effect of Promoter was the least of all significant factors in the model. The two MANOVA models showed similar results to those obtained with the mixed univariate models, except for iTarget x Genotype interaction. This interaction did not show a significant effect in the MANOVA analysis but it was significant in mixed models.

The NMDS ordination plot summarizes the results reported in this paper (Figure 2). The individuals are symbolized and colored according to the inverted repeat sequence used for gliadin silencing (‘g’, fragment targeted to γ-gliadins, ‘o’ chimeric fragment targeted to all three gliadins fraction, ‘go’, combination in different constructs of fragment ‘o’ and ‘g’). Individuals have been surrounded with ellipses to show the dispersion inside each iTarget. Besides, the centroids were plotted for each level of factor Year and Promoter. The most conspicuous feature of the plot was the clustering of individual samples according to the silencing fragment used. The ellipse that grouped together the control samples (iTarget C) was clearly separated from the rest. The same was observed for samples with the fragment targeted to γ-gliadin (iTarget g). Individuals that were transformed with the silencing fragment ‘o’ and the combination ‘go’, overlapped between themselves but none overlapped with control and ‘g’ silenced individuals. The ‘g’ silenced individuals had lower content of γ-gliadins, and a higher content of α-gliadins, ω-gliadins, total gliadins, LMW-GS, HMW-GS, total glutenins, and higher value of SDSS and Weight 1000 (thousand seed weight) than the controls lines. Individuals silenced with the ‘o’ fragment had a lower content of ω-gliadins, α-gliadins, total-gliadins, total prolamins and a lower ratio of Gli/Glu than the control and ‘g’ lines, though, like ‘g’ individuals, they had lower content of γ-gliadins too. Furthermore, these lines (‘o’) showed an increase of the content of LMW-GS, HMW-GS and total glutenin. Finally, individuals with silencing combination ‘go’ were very similar to those with only ‘o’, but the former had a lower overall glutenin content and a lower content of total prolamins. The Year and Promoter centroids did not show much variability, focusing around the center of the plot. However, the year 2010 showed a lower content of ω-gliadins, α-gliadins, total gliadins and total prolamins than year 2011. The Promoter factor levels showed little variation, with the Promoter ‘D’ being the most different from the other two levels (‘G’ y ‘GD’). The control samples and the silenced ‘g’ showed greater variability (greater dispersion of individuals samples) for content of ω-gliadins, α-gliadins and total gliadins.

thumbnailFigure 2. Non-metric multidimensional scaling (NMDS) graph depicting the relative importance of variables (arrows) explaining the transgenic and wild-type control lines composition. Ellipses denote the dispersion around the mean values of each iTarget. The individual points are the ordenation of each line and they are colored and symbolized according to their silencing level. The centroids are the average of the levels of the factors indicated. alpha, α-gliadins content; gamma, γ-gliadins content; omega, ω-gliadins content; gli, total gliadin content; HMW, total high molecular weight glutenin subunit content; LMW, total low molecular weight glutenin subunit content; glu, total glutenins content; gli.glu, ratio total gliadins content/total glutenins content; prolamin, prolamins content; non_gluten, non-gluten proteins; SDSml, SDS sedimentation test; Weight1000, thousand seeds weight.

Factor interactions analysis

The significant interactions in the mixed models reported above are shown graphically (Models 1 and 2) in Figure 3 and Additional file 1 for the variables with significant interaction between iTarget x Year and iTarget x Genotype.

thumbnailFigure 3. Interaction plots. Interaction of the iTarget with the genotype on (a) SDS sedimentation test (SDSS) values and (b) low molecular weight glutenin subunit content (LMW-GS). C, control lines; g, γ-gliadin silenced lines; o, ω/α-gliadin silenced lines; go, γ- and ω/α-gliadin silenced lines.

Additional file 1. Interaction plots of the iTarget with the year (environment) and genotype on the storage proteins contents and ratios. Only the significant interactions had been plotted. C, control lines; g, γ-gliadin silenced lines; o, ω/α-gliadin silenced lines; go, γ- and ω/α-gliadin silenced lines.

Format: PDF Size: 12KB Download file

This file can be viewed with: Adobe Acrobat ReaderOpen Data

Within iTarget x Genotype interactions, some of the most interesting were those for the SDSS test (Figure 3a) and the LMW-GS content (Figure 3b). In case of the SDSS test, the silenced lines with ‘g’ fragment of the genotype BW2003 had higher SDSS values than control lines but not those of genotype BW208, which showed equal SDSS values than control lines. For both genotypes, transgenic lines with combination of constructs ‘go’ had a marked decrease for the SDSS values in comparison to that of controls ‘C’ and ‘g’ lines. However, this decrease in the SDSS values was not as pronounced for genotype BW208 when the ‘o’ single fragment was used.

iTarget x Genotype interaction on the content of LMW-GS was decreased in lines with silencing ‘o’ of genotype BW2003 but increased in genotype BW208. It is important to highlight that the iTarget x Genotype interaction for the content of LMW-GS had the same pattern as that of SDSS, Total Prolamin and Total Glutenin (Figure 3a, Additional file 1k and i). However, iTarget x Genotype interaction on the content of HMW-GS was increased in both genotypes and iTarget (Additional file 1j).

As shown, for transgenic lines, year 2010 showed higher contents of ω-gliadins, α-gliadins, total gliadins and total prolamins but not of non-gluten proteins. This year effect was more prominent for the lines with silenced ‘g’ (Additional file 1a-f).

Discussion

In this paper we report 38 transgenic lines from two genotypes with different gliadin fractions silenced by RNAi, using three different promoter combinations and three different silencing construct combinations. The lines were tested during two consecutive years period to evaluate the effect of the promoters, silencing fragments and the environment on the down regulation of wheat gliadins and the effects of that silencing on other storage proteins and quality grain constituents.

Down-regulation of α-, γ-, and ω-gliadins by RNAi is an effective approach to reduce the expression of CD-related gliadin T-cells epitopes [15], which could be the basis for the development of products suitable not only for CD patients but also for other gluten intolerance patients. However, as consequence of this silencing, there is a re-balance of prolamin distribution, leading to the increment of total protein content in some particular lines [15,18], but not as a general effect. In most of the research articles where storage proteins of cereals were silenced by RNAi or mutation, the authors reported either a decrease or no variation in the total protein content [19-24] reported an overall increase of total protein in lines transformed with the ‘g’ antisense fragment, which agree with the data reported in this work. Although the total protein also increases in lines transformed with the silencing fragment ‘go’, this should be interpreted with caution, because the increase in total protein was associated with a decrease in thousand seed weight. It is known that grains that are not completely filled have higher protein content, and a higher embryo to endosperm protein ratio. Thus the increase in total protein in ‘go’ lines may come not only from non-gluten proteins (albumins and globulins) but also from a higher proportion of proteins from the embryo.

Lines transformed with silencing fragment ‘g’ showed a strong reduction of γ-gliadins, which was over-compensated by an increase of α-gliadins, ω-gliadins, HMW-GS and LMW-GS [24]. The term over-compensation is appropriate because the reduction of γ-gliadins results in a higher content of total protein rather than retaining the level of control lines. The reduction of all gliadins, with the silencing fragment ‘o’ or with the combination ‘go’, was accompanied by an increase in HMW-GS and a reduction of LMW’GS, except in line BW2003 with silencing ‘o’ where there was also an increase of LMW-GS. Increased HMW-GS was not enough to offset the lack of gliadins, resulting in the reduction of total content of prolamins, although the total protein remained constant or even increased in the case of the lines silenced with ‘go’ fragments. Therefore, the over-compensation of total protein must be from non-gluten proteins, such as albumins and globulins, as described previously [18] (and as shown by the calculation of non-gluten protein presented in this paper). This compensation effect had been reported in maize and rice, where the silencing of a group of storage proteins leads to a rearrangement of other storage proteins. In maize, the reduction in 22-kDa α-zeins levels by RNAi were compensated by increases of the 19-kDa α-zeins, and vice versa [25]. In the rice mutant line Low Glutelin Content-1 (LGC-1), the content of glutelin was reduced and the contents of other seed storage proteins, including prolamins, were increased [26]. Such up-regulation is not specific to the LGC-1 mutant and it is thought to be a non-specific compensation for the reduction of glutelin. On the other hand, reductions of glutelins and sulfur-rich 10-kDa prolamin levels by RNAi in rice were preferentially compensated by increases of sulfur-poor and other sulfur-rich prolamins, respectively, indicating that sulfur-containing amino acids might be involved in regulating seeds storage protein composition. It could be suggested that transgenic lines, with storage proteins reduced, attempt to compensate total protein content, first with related proteins, and then with unrelated proteins, if necessary. Therefore, the protein compensation could be a selective process because it does not use any kind of protein to compensate. If the compensation is governed by the availability of amino acids, the compensation process may be selectively determined by similarities in the amino acid composition of proteins.

The data presented in this paper allow the detection of differences between the two promoters used for gliadin silencing. Although the Promoter factor explained only part of the variability, it was clear that the γ-gliadin promoter had a higher efficiency which was demonstrated with a better silencing of γ-gliadins, but the efficiency was further increased when both promoters were used in combination. Both the γ-gliadin promoter and, the combination of promoters, led to a decrease of the content of LMW-GS. The contribution to the effectiveness of the promoters that drive a silencing fragment is determined by their expression level during the target expression [27,28]. Results reported by [29,30] concluded that D-hordein and γ-gliadin promoters both had high expression levels in the wheat endosperm but with different expression profiles. The D-hordein promoter was expressed in later stages of grain development than the γ-gliadin promoter. The high efficiency of the γ-gliadin promoter may be due to a higher expression level and/or to a better adjustment with the target genes. However, it is clear that not only the expression level is important, but also a broad expression profile which allows the expression of target genes throughout grain development. This may be what occurred when the combination of the two complementary promoters was used, leading to greater effectiveness.

Among the gliadin fractions, the content of ω-gliadins showed the highest variability being strongly environment-dependent (the factor year explained 69% of the variance). In fact it has been reported that the ω-gliadins modify their expression in response to S deficiency (increasing their expression with low S) [31-35], and N fertilization (the higher input of N the greater accumulation of ω-gliadins) [36-38]. Moreover, the proportions of ω-gliadins increase when grain is exposed to high temperature during grain filling [36,39]. In addition, the high environment-dependence of the ω-gliadin content also could be due to the fact that it is the group which is less efficiently silenced. Therefore, the iTarget factor has less influence on the ω-gliadin content compared with environmental factors.

The iTarget x Genotype interaction for SDSS showed that there was an increase of SDSS value in transgenic lines with iTarget ‘g’ but only for the BW2003 genotype, and a more pronounced decrease of SDSS for genotype BW2003 in comparison with genotype BW208 in iTarget lines ‘o’. Although the iTarget x Genotype interaction for SDSS in iTarget lines ‘o’ has first been analyzed in this work, the different behavior of genotypes with respect to SDSS has already been reported [24]. The latter interaction is associated in turn with the interaction of Genotype x iTarget for LMW-GS. In fact, in the BW2003 lines with a decrease of the SDSS values, this was associated with a decrease of LMW-GS. Moreover, the NMDS ordination graph showed that the LMW-GS content and SDSS value were highly correlated. A similar association between quality parameters and LMW-GS content was also reported by [24], who showed a positive correlation between some mixograph parameters and the SDSS test with an individual LMW-GS peak and total LMW-GS contents.

Conclusions

iTarget was the main factor that affects the characteristics measured in the the transgenic and control lines analyzed. The silencing is a stable effect over plant generations considered in this study, and produces an efficient and lasting reduction of the different gliadin fractions. Although, the gliadin silencing efficiency is determinated mainly by the factor iTarget, but the promoter, the genotype and the environment factors also affect the gliadin silencing and these factors should be taken into account to achieve a gliadin silencing as high as possible.

The ω-gliadins are the prolamins with the greatest variability, and they are less efficiently silenced. Future work should be directed to the manipulation of this protein fraction to obtain greater and more robust reductions in its content. γ-gliadin silencing is compensated with an increase of α- and ω-gliadins, and to a lesser extent with glutenins. The silencing of all gliadins is compensated by an increase of the glutenin content but also with an increase of non-gluten proteins.

The different promoters tested used alone have a similar efficiency, but the combination of the two promoters with different patterns of expression allows us to cover a wider range of grain developmental stages and thereby higher silencing efficiency. This, may be a general strategy for efficient gene silencing of gene families like gliadins in which constructs under two promoters with different expression profiles, covering as much as possible the expression of target genes during all stages of grain filling.

The results presented will underpin the development of future strategies to achieve more precise and effective silencing of toxic peptides in relation to reducing gluten intolerance pathologies.

Methods

Plant material

Twenty-two transgenic lines of T. aestivum cv ‘Bobwhite 208’ (BW208) and sixteen transgenic lines of T. aestivum cv ‘Bobwhite 2003’ (BW2003) and their corresponding wild-type lines were used in this study. Line BW2003 carries the translocation T1BL.1RS from rye. Four hairpin RNA (hpRNA) vectors were used to down regulate the γ-, α- and ω-gliadins: the pghpg8.1 and pGghpg8.1 vectors down-regulate the γ-gliadins; pDhp- ω/α and pGhp- ω/α vectors down-regulate all gliadin fractions [15]. The constructs pghpg8.1 and pDhp- ω/α contain the D-hordein promoter whereas the constructs pGghpg8.1 and pGhp- ω/α contain the γ-gliadin promoter [29,30]. Transgenic lines A1152, A1158, A1406, C655, C657, D445, D623, C217 and D598 contain the pghpg8.1 vector; lines D577, D682, D715, D716, D815, 22A, 22C, 24A and 24C contain the pGhpg8.1 vector; lines 28A, 28B, D783, E140 and E146 contain the pDhp- ω/α vector; lines D770, D793 and D894 contain the pGhp- ω/α vector; lines D874 and D876 contain pDhp- ω/α and pGhp- ω/α vectors; lines E33, E35, E39, E42, E76, E82, E83 and E122 contain both the pghpg8.1 and pDhp- ω/α vector; lines E93 and E96 contain the pGghpg8.1 and pGhp- ω/α vectors (Table 3). All transgenic lines were previously reported or obtained as described by [15,24] and self-pollinated for four generations to obtain homozygous lines. The seeds used in the assay of second year (year 2011) were from self-pollinated plants of the first assay (year 2010).

Reversed-phase high-performance liquid chromatography (RP-HPLC)

Gliadins and glutenins were extracted and quantified by RP-HPLC following the protocol reported by [24].

Thousand seed weight, total protein and starch, sodium dodecyl sulphate sedimentation (SDSS) assay, and non-gluten proteins content calculation

Thousand seed weight (g) was determined for 1000 seeds from each sample. Two measurements were carried out for each sample.

The protein content of whole flour was calculated from the Kjeldahl nitrogen content (%N × 5.7) according to the standard ICC method no. 105/2 [40]. Starch content was determined according to the standard ICC method no. 123/1 [41]. Both parameters were expressed on a dry matter basis.

The SDS sedimentation volume was determined as described by [42]. Three replicates were carried out for each experimental unit.

The non-gluten proteins, expressed in μg/mg, were calculated as follow: [Total protein (%) * ((100 - moisture in%)/100) - (Prolamins content in μg/mg/10)] * 10.

Experimental design and statistical analysis

All analyses and plots were conducted with the statistical software R version 2.14.1 [43]. The experimental design was a RCBD with three replications of each line and five plants per plot. This experimental design was repeated two consecutive years (2010 and 2011). The randomized block designs were generated with the package agricolae [44]. The variables were submitted to a linear mixed model fitted by restricted maximum likelihood (REML) using function lmer [45]. The models were adjusted using the factors ‘Year’, ‘Block’, ‘Genotype’, ‘iTarget’, ‘ Promoter’ and ‘Line’ and the covariate ‘Total protein’, where ‘Year’ is the year in which the assays were performed (2010 and 2011), ‘Block’ are the three blocks of each RCBD, ‘Genotype’ is the genotype of transgenic and wild types lines (BW208 and BW2003), ‘iTarget’ is the silencing fragment used to silent the gliadins (‘C’ is the control, ‘g’ is the construct to silence γ-gliadins, ‘o’ is the construct to silence all the gliadins, ‘go’ is when the two constructs are used in combination), ‘Promoter’ is the promoter used to drive the silencing fragments (‘D’ is the D-hordein promoter and ‘G’ is the γ-gliadin promoter), and ‘Line’ represent the different transgenic and controls lines used. Factors ‘Year’ and ‘Block’ were considered of random effects, whereas ‘iTarget’, ‘Promoter’ and ‘Line’ were considered of fixed effects. ‘Total protein’ was used as fixed covariate because grain protein content is environment dependent and it has a strong effect on the storage proteins. The lmer model used for testing the effect of ‘Year’, ‘Total protein’, ‘Genotype’ and ‘iTarget’ was: ‘variable ~ Total protein + (1|Year) + (1|Block:Year) + (1|iTarget:Year) + (1|Genotype:Year) + iTarget*Genotype’ (Model 1). Data were normalized against their respective controls before testing the effect of the two promoters (‘D’ and ‘G’) to drive the expression of the silencing fragments. To check the effect of the factor ‘Promoter’, the following lmer model with the data normalized was used: ‘variable ~ Total protein + (1|Year) + (1|Block:Year) + (1|iTarget:Year) + (1|Promoter:Year) + iTarget + Promoter’ (Model 2). The differences between the control and the transgenic lines were assessed using the following mixed effect model: ‘variable ~ Total protein + (1|Year) + (1|Block:Year) + Line’ (Model 3). Residuals were tested for normal distribution and for homogeneity of variances using the model criticism plots generated by the function mcp.fnc (package LMERConvenienceFunctions) [46]. In the cases where the conditions of normality and homogeneity of variances were violated, the Box-Cox transformation was applied (function powerTransform, package car; [47]). p-values for the analysis of variance (or deviance) as well as the amount of deviance explained (%) for each fixed-effect of mixed models were calculated with the function pamer.fnc (package LMERConvenienceFunctions). Post hoc multiple-comparison was carried out with the function glht (package multcomp; [48]).

To identify the relative importance of each factor in the total set of data, we performed a non-parametric multivariate analysis of variance (MANOVA) using the function adonis (package Vegan; [49]). The non-parametric MANOVA model was: ‘variable ~ Total protein + Block%in%Year + Year*Genotype + Year*iTarget + Genotype*iTarget + Line’ (Model M1). The effect of the two promoters in the data set was compared using the following model with normalized data: ‘variable ~ Total protein + Block%in%Year + Year*iTarget + Year*Promoter + Promoter*iTarget + Line’ (Model M2). Nonmetric Multidimensional Scaling (NMDS) ordination plots were generated to demonstrate trends in transgenic and control lines in relation to the variables studied using the function metaMDS (package Vegan). For both analyses (adonis and metaMDS), the method used to calculate the pairwise distances was ‘Euclidean’ with 999 permutations. The package maptools was used to avoid the label overlapping in the ordination plot [50].

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

FB conceived of the study, designed experiments, interpreted the results and wrote the paper. JG conceived of the study, designed experiments, interpreted the results and wrote the paper. FP conceived of the study, designed and carried out experiments, analyzed and interpreted the results, and wrote the paper. All authors read and approved of the final manuscript.

Acknowledgments

The Spanish Ministry of Science and Innovation (AGL2010-19643-C02-02 and TRA2009_0047), the European Regional Development Fund (FEDER), and the Junta de Andalucía (Project P09AGR-4783) supported this work. Fernando Pistón thanks the Spanish program Ramón y Cajal. Javier Gil-Humanes acknowledges financial support from the I3P Program from the Agencia Estatal Consejo Superior de Investigaciones Científicas, which is co-financed by the European Social Fund. We thank Dr. Paul Lazzeri (Agrasys SL) for critical comments on the manuscript and English editing. The technical assistance of Ana García is also acknowledged.

References

  1. Shewry PR, Miles MJ, Tatham AS: The prolamin storage proteins of wheat and related cereals.

    Prog Biophys Mol Biol 1994, 61:37-59. PubMed Abstract OpenURL

  2. Shewry PR, Halford NG: Cereal seed storage proteins: structures, properties and role in grain utilization.

    J Exp Bot 2002, 53:947-958. PubMed Abstract | Publisher Full Text OpenURL

  3. Zuidmeer L, Goldhahn K, Rona RJ, Gislason D, Madsen C, Summers C, Sodergren E, Dahlstrom J, Lindner T, Sigurdardottir ST, McBride D, Keil T: The prevalence of plant food allergies: a systematic review.

    J Allergy Clin Immunol 2008, 121:1210-1218.

    e4

    PubMed Abstract | Publisher Full Text OpenURL

  4. Mariné M, Farre C, Alsina M, Vilar P, Cortijo M, Salas A, Fernández-Bañares F, Rosinach M, Santaolalla R, Loras C, Marquès T, Cusí V, Hernández MI, Carrasco A, Ribes J, Viver JM, Esteve M: The prevalence of coeliac disease is significantly higher in children compared with adults.

    Aliment Pharmacol Ther 2011, 33:477-486. PubMed Abstract | Publisher Full Text OpenURL

  5. Sapone A, Lammers KM, Casolaro V, Cammarota M, Giuliano MT, Rosa MD, Stefanile R, Mazzarella G, Tolone C, Russo MI, Esposito P, Ferraraccio F, Cartenì M, Riegler G, De ML, Fasano A: Divergence of gut permeability and mucosal immune gene expression in two gluten-associated conditions: celiac disease and gluten sensitivity.

    BMC Med 2011, 9:23. PubMed Abstract | BioMed Central Full Text | PubMed Central Full Text OpenURL

  6. Sollid LM: Coeliac disease: dissecting a complex inflammatory disorder.

    Nat Rev Immunol 2002, 2:647-655. PubMed Abstract | Publisher Full Text OpenURL

  7. Molberg Ø, Solheim Flaete N, Jensen T, Lundin KEA, Arentz-Hansen H, Anderson OD, Kjersti Uhlen A, Sollid LM: Intestinal T-cell responses to high-molecular-weight glutenins in celiac disease.

    Gastroenterology 2003, 125:337-344. PubMed Abstract | Publisher Full Text OpenURL

  8. Arentz-Hansen H, Mcadam SN, Molberg Ø, Fleckenstein B, Lundin KEA, Jørgensen TJD, Jung G, Roepstorff P, Sollid LM: Celiac lesion T cells recognize epitopes that cluster in regions of gliadins rich in proline residues.

    Gastroenterology 2002, 123:803-809. PubMed Abstract | Publisher Full Text OpenURL

  9. Sollid LM, McAdam SN, Arentz-Hansen H, Korner R, Molberg O, Quarsten H, Vader W, Kooy YMC, Lundin KEA, Koning F, Roepstorff P: The intestinal T cell response to α-gliadin in adult celiac disease is focused on a single deamidated glutamine targeted by tissue transglutaminase.

    J Exp Med 2000, 191:603-612. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  10. Payne PI: Genetics of wheat storage proteins and the effect of allelic variation on bread-making quality.

    Ann Rev Plant Physio 1987, 38:141-153. OpenURL

  11. Anderson OD, Greene FC: The α-gliadin gene family. II. DNA and protein sequence variation, subfamily structure, and origins of pseudogenes.

    Theor Appl Genet 1997, 95:59-65. Publisher Full Text OpenURL

  12. Sabelli PA, Shewry PR: Characterization and organization of gene families at the Gli-1 loci of bread and durum wheats by restriction fragment analysis.

    Theor Appl Genet 1991, 83:209-216. OpenURL

  13. Baulcombe D: RNA silencing in plants.

    Nature 2004, 431:356-363. PubMed Abstract | Publisher Full Text OpenURL

  14. Gil-Humanes J, Pistón F, Hernando A, Alvarez JB, Shewry PR, Barro F: Silencing of γ-gliadins by RNA interference (RNAi) in bread wheat.

    J Cereal Sci 2008, 48:565. Publisher Full Text OpenURL

  15. Gil-Humanes J, Pistón F, Tollefsen S, Sollid LM, Barro F: Effective shutdown in the expression of celiac disease-related wheat gliadin T-cell epitopes by RNA interference.

    Proc Natl Acad Sci USA 2010, 107:17023-17028. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  16. Burnouf T, Bietz JA: Chromosomal control of glutenin subunits in aneuploid lines of wheat: analysis by reversed-phase high-performance liquid chromatography.

    Theor Appl Genet 1985, 70:610-619. Publisher Full Text OpenURL

  17. Wieser H, Antes S, Seilmeier W: Quantitative determination of gluten protein types in wheat flour by reversed-phase high-performance liquid chromatography.

    Cereal Chem 1998, 75:644-650. Publisher Full Text OpenURL

  18. Gil-Humanes J, Pistón F, Shewry PR, Tosi P, Barro F: Suppression of gliadins results in altered protein body morphology in wheat.

    J Exp Bot 2011, 62:4203-4213. PubMed Abstract | Publisher Full Text OpenURL

  19. Hartings H, Lauria M, Lazzaroni N, Pirona R, Motto M: The Zea Mays mutants opaque-2 and opaque-7 disclose extensive changes in endosperm metabolism as revealed by protein, amino acid, and transcriptome-wide analyses.

    BMC Genomics 2011, 12:41. PubMed Abstract | BioMed Central Full Text | PubMed Central Full Text OpenURL

  20. Wu Y, Messing J: RNA interference can rebalance the nitrogen sink of maize seeds without losing hard endosperm.

    PLoS ONE 2012, 7:e32850. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  21. Schmidt MA, Barbazuk WB, Sandford M, May G, Song Z, Zhou W, Nikolau BJ, Herman EM: Silencing of soybean seed storage proteins results in a rebalanced protein composition preserving seed protein content without major collateral changes in the metabolome and transcriptome.

    Plant Physiol 2011, 156:330-345. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  22. Lange M, Vincze E, Wieser H, Schjoerring JK, Holm PB: Suppression of C-hordein synthesis in barley by antisense constructs results in a more balanced amino acid composition.

    J Agr Food Chem 2007, 55:6074-6081. Publisher Full Text OpenURL

  23. Kawakatsu T, Hirose S, Yasuda H, Takaiwa F: Reducing rice seed storage protein accumulation leads to changes in nutrient quality and storage organelle formation.

    Plant Physiol 2010, 154:1842-1854. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  24. Pistón F, Gil-Humanes J, Rodríguez-Quijano M, Barro F: Down-regulating γ-gliadins in bread wheat leads to non-specific increases in other gluten proteins and has no major effect on dough gluten strength.

    PLoS ONE 2011, 6:e24754. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  25. Huang S, Adams WR, Zhou Q, Malloy KP, Voyles DA, Anthony J, Kriz AL, Luethy MH: Improving nutritional quality of maize proteins by expressing sense and antisense zein genes.

    J Agr Food Chem 2004, 52:1958-1964. Publisher Full Text OpenURL

  26. Iida S, Amano E, Nishio T: A rice (oryza sativa L.) mutant having a low content of glutelin and a high content of prolamine.

    Theor Appl Genet 1993, 87:374-378. Publisher Full Text OpenURL

  27. Kerschen A, Napoli CA, Jorgensen RA, Müller AE: Effectiveness of RNA interference in transgenic plants.

    FEBS Lett 2004, 566:223-228. PubMed Abstract | Publisher Full Text OpenURL

  28. Marjanac G, Karimi M, Naudts M, Beeckman T, Depicker A, De Buck S: Gene silencing induced by hairpin or inverted repeated sense transgenes varies among promoters and cell types.

    New Phytol 2009, 184:851-864. PubMed Abstract | Publisher Full Text OpenURL

  29. Pistón F, Shewry PR, Barro F: D hordeins of hordeum chilense: a novel source of variation for improvement of wheat.

    Theor Appl Genet 2007, 115:77-86. PubMed Abstract | Publisher Full Text OpenURL

  30. Pistón F, Marín S, Hernando A, Barro F: Analysis of the activity of a γ-gliadin promoter in transgenic wheat and characterization of gliadin synthesis in wheat by MALDI-TOF during grain development.

    Mol Breeding 2009, 23:655-667. Publisher Full Text OpenURL

  31. MacRitchie F, Gupta R: Functionality-composition relationships of wheat flour as a result of variation in sulfur availability.

    Aust J Agr Res 1993, 44:1767-1774. Publisher Full Text OpenURL

  32. Moss H, Wrigley C, MacRichie R, Randall P: Sulfur and nitrogen fertilizer effects on wheat. II. Influence on grain quality.

    Aust J Agr Res 1981, 32:213-226. Publisher Full Text OpenURL

  33. Wieser H, Gutser R, Von Tucher S: Influence of sulphur fertilisation on quantities and proportions of gluten protein types in wheat flour.

    J Cereal Sci 2004, 40:239-244. Publisher Full Text OpenURL

  34. Wrigley CW, Du Cros DL, Fullington JG, Kasarda DD: Changes in polypeptide composition and grain quality due to sulfur deficiency in wheat.

    J Cereal Sci 1984, 2:15-24. Publisher Full Text OpenURL

  35. Zhao F-J, Salmon SE, Withers PJA, Evans EJ, Monaghan JM, Shewry PR, McGrath SP: Responses of breadmaking quality to sulphur in three wheat varieties.

    J Sci Food Agr 1999, 79:1865-1874. Publisher Full Text OpenURL

  36. Daniel C, Triboi E: Effects of temperature and nitrogen nutrition on the grain composition of winter wheat: effects on gliadin content and composition.

    J Cereal Sci 2000, 32:45-56. Publisher Full Text OpenURL

  37. Dupont FM, Altenbach SB: Molecular and biochemical impacts of environmental factors on wheat grain development and protein synthesis.

    J Cereal Sci 2003, 38:133-146. Publisher Full Text OpenURL

  38. Wieser H, Seilmeier W: The influence of nitrogen fertilisation on quantities and proportions of different protein types in wheat flour.

    J Sci Food Agr 1998, 76:49-55. Publisher Full Text OpenURL

  39. Corbellini M, Canevar MG, Mazza L, Ciaffi M, Lafiandra D, Borghi B: Effect of the duration and intensity of heat shock during grain filling on dry matter and protein accumulation, technological quality and protein composition in bread and durum wheat.

    Funct Plant Biol 1997, 24:245-260. OpenURL

  40. ICC: Determination of crude protein in cereals and cereal products for food and for feed.

    Int Assoc Cereal Sci Technol Method 1994, 105:2. OpenURL

  41. ICC: Determination of starch content by hydrochloric acid dissolution.

    Int Assoc Cereal Sci Technol Method 1994, 123:1. OpenURL

  42. Williams P, El-Haramein FJ, Nakkoul H, Rihawi S: Crop quality evaluation methods and guidelines. Syria: ICARDA Technical Manual, International Center for Agricultural Research in the Dry Areas; 1988. OpenURL

  43. Team RDC: R: a language and environment for statistical computing. Vienna, Austria; 2012. OpenURL

  44. de Mendiburu F:

    Agricolae: statistical procedures for agricultural research. 2012. OpenURL

  45. Bates D, Maechler M, Bolker B:

    lme4: linear mixed-effects models using S4 classes. 2011. OpenURL

  46. Tremblay A:

    LMERConvenienceFunctions: a suite of functions to back-fit fixed effects and forward-fit random effects, as well as other miscellaneous functions. 2012. OpenURL

  47. Fox J, Weisberg S: An R companion to applied regression. 2nd edition. Thousand Oaks CA: Sage; 2011. OpenURL

  48. Hothorn T, Bretz F, Westfall P: Simultaneous inference in general parametric models.

    Biometrical J 2008, 50:346-363. Publisher Full Text OpenURL

  49. Oksanen J, Blanchet FG, Kindt R, Legendre P, Minchin PR, O’Hara RB, Simpson GL, Solymos P, Stevens MHH, Wagner H: Community ecology package. Vegan; 2012. OpenURL

  50. Lewin-Koh NJ, Bivand R, EJ P c b, Archer E, Baddeley A, Bibiko H-J, Callahan J, Dray S, Forrest D, Friendly M, Giraudoux P, Golicher D, Rubio VG, Hausmann P, Hufthammer KO, Jagger T, Luque SP, MacQueen D, Niccolai A, Short T, Snow G, Stabler B, Turner R:

    Maptools: tools for reading and handling spatial objects. 2012. OpenURL