Arabidopsis thaliana accessions are naturally occurring and essentially homozygous inbred lines that are frequently used to investigate genetic and/or metabolic variations [1–4] and to identify genes relevant for intraspecific adaptation phenomena in plants [5–7]. Due to intraspecific genetic variation, many A. thaliana accessions differ in growth and development even when grown alongside under the same conditions [8, 9]. Since their divergence from A. lyrata 5-10 million years ago, A. thaliana accessions possess a long evolutionary history of intraspecific diversification [10, 11]. Genetic variation leads to nucleotide polymorphisms in both coding and noncoding gene regions. Nonsynonymous substitutions locally alter the amino acid sequence of either the transit peptide or the mature protein at the level of translation and, thereby, potentially may affect protein-related functions. Synonymous substitutions do not alter the amino acid sequence but may affect level and/or stability of the transcripts as well as the rate of translation. Thereby, they might indirectly alter the level of a given protein. Likewise, variation in noncoding regions is not translated into amino acid polymorphisms, but can exert diverse effects, such as alternative splicing, introduction of premature stop codons of transcription or translation, altered transcripts stability and/or rate of gene expression [12]. Thereby, it indirectly may also affect level and/or amino acid sequence of a given protein.

To a large extent, regulation of gene expression is based on the action of regulatory elements that are located in positions designated as cis (i.e. close to the target gene) or trans (i.e. distant from the gene) [13, 14]. In summary, genetic diversity may affect phenotypic traits by acting on different levels ranging from gene expression to transcript level and altered features and/or functions of the protein.

In many A. thaliana accessions studied so far, vegetative biomass is negatively correlated with leaf starch content [15]. Therefore, intraspecific genetic variation appears to massively affect the central carbon metabolism and growth of the entire plant. In A. thaliana, as in many other plant species, transitory starch is a major product of photosynthesis which is deposited in the stromal space of the mesophyll chloroplasts as water-insoluble particles (i.e. granules). They possess a strictly defined and evolutionary conserved (inter)molecular order [16, 17] and consist of two types of polyglucans, amylopectin and amylose. Amylopectin is a large, highly branched polyglucan representing the main constituent of the granule. By contrast, amylose is a polydispers and essentially unbranched biopolymer that, in most cases, is a minor starch compound and contributes little to the (inter)molecular organization of the entire starch particle [17].

Starch biochemistry is based on the coordinated and evolutionary conserved action of 30 to 40 (iso)enzymes and is more complex than the classical glycogen metabolism [18]. Massive starch biosynthesis proceeds by a repetitive glucosyl transfer from an appropriate donor (such as ADPglucose) to non-reducing ends of oligo- or polyglucans [19]. ADPglucose-dependent chain elongation is mediated by at least five classes of starch synthases (ADP-Glc: α-1,4 glucan α-4-glucosyl transferase; EC 2.4.1.21). In A. thaliana, each class is represented by only a single gene. Based on sequence similarity, kinetic properties, and the occurrence of consensus motifs, they are all related to the glycogen synthases from both prokaryotes and eukaryotes and are members of the glucosyl transferase family 5 (GT5).

Unlike in A. thaliana, in some lower and higher plants starch synthase classes are represented by more than a single gene. Complexity of the starch synthase classes tends to increase if cells or tissues are capable of metabolizing various starch pools (such as transitory and reserve starches) that are spatially and/or temporarily separated (Table 1).

The five SS classes (and the glycogen synthases as well) share a core region of approximately 60 kDa that is indispensable for catalytic activity and consists of two conserved domains often designated as GT5 (glycosyl transferase family 5) and GT1 (glycosyl transferase family 1). Both domains are separated by a short and more variable linker region. The GT5 domain is typical for glucosyl transferases following a retaining mechanism [25]. GT1 is consistently located close to the C-terminus of the starch/glycogen synthases and is frequently found in glycosyl transferases mediating an inverting mode of glycosyl transfer [26]. In all starch synthases, both the GT5 and GT1 domain are involved in binding of the glucosyl donor, ADPglucose, and together they form the catalytically active region of starch synthases.

A scheme of the domain structure of the five A. thaliana starch synthases and of three prokaryotic glycogen synthases is presented in Figure 1. AtGBSS largely consists of the core region containing the two domains GT5 and GT1 and possesses only a short additional sequence at the C-terminus. The size of AtGBSS is similar to that of the prokaryotic glycogen synthases. All soluble starch synthases (AtSSI to AtSSIV) carry a N-terminal extension that might modulate catalytic activity. The length of the extension strongly varies both among SS classes of a given plant species and between species. A unique feature of the N-terminal sequence of AtSSII is a serine-rich region of uncertain function. AtSSIII and AtSSIV possess large N-terminal extensions that are similar in size, but diverse in amino acid sequence. The N-terminal sequence of AtSSIII contains three repeats of a distinct carbohydrate binding module (CBM25) [27, 28] that are not present in AtSSIV (Figure 1).

In previous studies, AtSSI has been localized in one of the few genome regions that possess a high intraspecific variability [29]. Here, we used 30 A. thaliana accessions and sequenced genes encoding the five starch synthases including their promoter regions. Accessions were selected on the basis of both largely different climates at their original habitats and wide geographic distribution of the natural populations. The objective was to estimate the gene-specific level of variation as well as possible implications for gene expression and protein structure. Specifically, we were interested in: (i) whether genes of the five starch synthase classes exhibit a similar degree of both synonymous and nonsynonymous variation, (ii) whether specific gene trees of the five starch synthase classes show the same pattern of haplotype clustering across accessions, (iii) whether selection can be inferred to act on genes and/or single polymorphic sites, (iv) whether polymorphisms across accessions and starch synthase classes have functional implications, and (v) whether there is a relationship between genomic variation and transcript levels in starch synthases.

Primary metabolism has been defined as ‘those essential reactions involving compounds that are formed as part of the normal anabolic and catabolic processes, which result in assimilation, respiration, transport, and differentiation processes that take place in most, if not all, cells of an organism’ [39]. It has been assumed that genes involved in primary metabolism are more conserved than secondary metabolism genes because of their essential function [40]. For this reason, gene-specific investigations of intraspecific variability among accessions of A. thaliana were largely performed with secondary metabolism genes [41–49]. Currently, numerous whole genomes become available as part of the 1,001 genomes project [50] become available. However, with these whole-genome approaches subtle differences within and across highly homologous gene loci can be overlooked and are easier detectable when targeting specific genes and their adjacent genomic regions [51]. Here, we demonstrate that genes encoding the starch synthases and exert essential functions in the plant primary metabolism possess high levels of nucleotide diversity as genes related to the secondary metabolism [41–49]. Furthermore, we show that transcript levels of starch synthases vary among accessions in all starch synthases, but both amount and pattern of variation differ between starch synthases. Variation was minor in AtSSII, AtSSIII, and AtGBSS, but higher in AtSSI, and AtSSIV. The differences in transcript levels clearly correlate with the gene tree in AtSSIV. We argue that such a strict correlation between combined promoter/gene sequences and transcript levels is indicative of cis-regulation. Such correlation was absence in AtSSI, AtSSII, AtSSIII, and AtGBSS, indicating trans-regulation is the major regulatory mechanism in these genes.

Obviously, functional analyses are expected to be more complete if the genetic variations were compared with the total starch synthase activity and/or the zymograms obtained for the various accessions. Unfortunately, this approach is not possible. Some starch synthases apparently contribute very little to the total enzyme activity measured in crude extracts. Furthermore, in zymograms performed with leaf extracts some AtSS isozymes are recovered as multiple bands but products of other AtSS genes are not detectable at all (although the respective recombinant proteins exhibit enzyme activity). Thus, zymograms do not reflect the genetic complexity of the starch synthases [52].

With regard to the specific pattern of variation, our study is able to address the following questions:

(i) Do genes of the five starch synthase classes exhibit a similar degree of both synonymous and nonsynonymous variation?

The five starch synthases exhibit different degrees of both synonymous and nonsynonymous variation. AtSSI, AtSSIII, and AtGBSS are the most variable genes exhibiting similar nucleotide diversities. AtSSII and AtSSIV possess the lowest nucleotide diversity. We confirm that AtSSI is highly diverse [29]. Presumably, nucleotide diversity correlates with the in vivo function of the protein. As outlined above, the five starch synthase classes are likely to exert non-identical in vivo functions. SSI to SSIII classes are involved in amylopectin biosynthesis, whereas GBSS is essential for the biosynthesis of amylose. AtSSIV appears to be essential in one route of the initiation of starch granule biosynthesis [30, 31]. Mutants from A. thaliana lacking functional AtSSIV possess only a single, enlarged granule per plastid of the mesophyll cells [31]. Furthermore, overexpression of AtSSIV leads to increasing levels of both transitory and storage starch [32]. This specific function might be the reason for the lower nucleotide diversity in AtSSIV. Knockout mutants of AtSSII are deficient in intermediate α-glucan chain length with DP (degree of polymerization) of 12–25 [21]. This is compared to AtSSI (DP 8–12; [21]) and AtSSIII (DP 14–20; [53]) a broader spectrum. The lower level of genetic variation in AtSSII might be due to the partly overlapping function, by which AtSSII can substitute either AtSSI or AtSSIII. However, the reason for the relatively high level of nucleotide diversity in AtGBSS is not clear, as its ability to synthesize amylose is an important biological function [17].

(ii) Do specific gene trees of the five starch synthase classes show the same pattern of haplotype clustering across accessions?

One could assume that functionally (originated from an ancestral gene) related genes exhibit similar gene trees. In general, starch synthases show at least two haplogroups, except for AtSSII. The existence of diverse clusters of haplotypes has also been reported from Rubisco genes [51]. These divergent haplogroups are indicative of a relatively large long-term effective population size of the species and are likely to comprise ancient standing variation, rather than in situ divergence among accessions. As in Rubisco [51], the phylogenetic trees of the different starch synthase genes are neither congruent with one another nor do they reflect any geographical or ecological pattern. Because accessions are mainly homozygous inbred lines and, due to the local distance between naturally occurring populations, any exchange of gene variants between accessions is unlikely to occur. The incongruent phylogenetic pattern, however, could be indicative of relatively frequent recombination across accession (on an evolutionary timescale), by which gene tree and species/accession tree are disentangled. Only if many loci in combination or even whole genomes are phylogenetically analyzed, a reliable phylogeographic pattern can be detected among A. thaliana accessions [54, 55].

(iii) Can selection be inferred to act on genes and/or single polymorphic sites?

The Z-tests revealed purifying selection is acting on AtSSI, AtSSIII, AtSSIV, and AtGBSS. Previous studies on AtSSI yielded significant positive values for Tajima’s D, interpreted as indication for balancing selection [29, 56]. We also obtained positive (yet statistically not significant) D values for AtSSI. Because A. thaliana accessions comprise essentially homozygous inbred lines, we argue that such positive selection across accessions (as indicated by positive D values) should be interpreted as disruptive selection, i.e., divergent evolution among evolutionary lineages. As all the different accessions – because of their different geographic origin and a high level of selfing - very likely do not contribute to a single common gene pool, negative frequency-dependent (= balancing) selection across accessions (as postulated in [56]) appears less likely to us. Otherwise, a scenario of divergent evolution at AtSSI fits well the description of two haplotypes [56], a Col-0 type (designated as A type) and the Ler type (B type). According to our analysis, these groups are separated by 105 SNPs (thereof 21 in coding sequences) in AtSSI.

(iv) Do polymorphisms across accessions and starch synthase classes have functional implications?

Several nonsynonymous substitutions were found in each of the starch synthases, some of which located at positions that are highly conserved among plant species. For each starch synthase, several amino acids have been identified that are of particular functional importance [27, 28, 34, 57]. Furthermore, we performed an interspecies comparison and searched for polymorphisms at sites that are involved in ADP-glucose binding as well as catalysis by generating an alignment including protein sequences of all available starch synthases from maize, rice, and A. thaliana. None of those highly conserved sites were substituted in any of the accessions analyzed in this study. During starch biosynthesis many enzyme activities closely cooperate and, therefore, any disturbance of this concerted action may result in complex alterations of the starch structure. Several starch-related enzymes are likely to undergo protein-protein interactions in vivo and, therefore, may exert their biochemical functions mainly (or exclusively) as constituent of a protein complex rather than as a single catalyst [58, 59]. For heterotrophic tissues starch-related protein complexes have been described that consist of distinct starch synthases and branching enzyme (and the plastidial phosphorylase as well) and the formation of these high molecular weight complexes appears to be controlled by covalent protein modifications [58, 59]. If a given enzyme is inefficient or even not functional due to amino acid substitutions, the resulting functional implications may reflect not only that of a single enzyme, but rather that of the respective protein complex. The lack of polymorphism at functionally crucial sites and the inferred pattern of purifying selection hence underline the functional importance and evolutionary conservation of these genes. It is highly unlikely that new mutation will be positively selected, although such mutations could – in theory - establish partly or completely novel complexes or functionalities.

The promoters of starch synthases exhibit numerous polymorphisms and indels, which potentially influence the transcript level. We found significant differences in transcript levels among accessions for each starch synthase. In AtSSIV, the only starch synthase gene inferred to be cis-regulated (see above), we found several positions in the promoter region which show a diversity pattern congruent with haplogroup assignment and transcript level. However, these motifs could not be assigned to any known regulatory element in the PPDB.

(v) Is there a relationship between genomic variation and transcript levels in starch synthases?

Analyses of differences in transcript levels between accessions and their position in the combined promoter and gene phylogenetic tree revealed different kinds of relationships. In AtSSII, AtSSIII, and AtGBSS, we were unable to detect a tight correlation between genetic variation and transcript levels. We propose that trans-regulation is the major actor in these genes, because the polymorphisms that are responsible for differences in gene expression are apparently far away from the gene [13, 14]. By contrast, in AtSSIV, we found a clear correlation between genomic polymorphisms and transcript levels. For this reason, we propose cis-regulation as a major actor in AtSSIV expression [13, 14]. AtSSI is somehow peculiar. We could not detect any differences in transcript levels among accessions representing haplotype cluster I. It seems that none of the polymorphisms found in the promoters of these accessions occurred in functionally relevant cis-elements. In contrast, in accessions of cluster II cis-regulatory elements might be affected by polymorphisms, because we could detect expression variation among them. We assume that both cis- as well as trans-regulation are active in the regulation of the expression of AtSSI.

AtSSI, AtSSII, and AtSSIII are important for building the amylopectin chains. They possess partly overlapping, but also individual properties in α-glucan chain elongation, while GBSS is mainly responsible for amylose synthesis [17, 21, 30, 53, 60].

While screening 30 accessions of A. thaliana gene specifically, we detected several nonsynonymous substitutions in each of the five analyzed starch synthases (AtSSI, AtSSII, AtSSIII, AtSSIV, AtGBSS). Gene trees for single genes often revealed a clear-cut clustering of accessions, which is – however – not consistent across different starch synthase genes. Our results are compatible with previous findings that two haplogroups might adaptively (by positive selection) diverge in AtSSI. In general, we identified amino acid substitutions in the catalytic glucosyl transferase domains (GT5, GT1) in almost all enzymes. Some of these observed amino acid substitutions affect sites known to be highly conserved across different plant species. Transcript analyses revealed significant differences in all starch synthases, although the extent varies among them. Comparison of haplotype clustering and transcript levels of starch synthases is indicative of both trans- and cis-regulated genes. AtSSI, AtSSII, AtSSIII, and AtGBSS are important for chain elongation of amylopectin and amylose, respectively. The corresponding genes are mainly globally regulated by trans-regulation, whereby elongation of α-glucan chains might be coordinated in concert. AtSSIV has exclusive and far reaching function and a separate regulation by individual cis-regulatory transcription factors – as inferred by our analysis - appears reasonable.