β-tubulin paralogue tubC is frequently misidentified as the benA gene in Aspergillus section Nigri taxonomy: primer specificity testing and taxonomic consequences

β-tubulin (benA, tub-2) and calmodulin (caM) are crucial genes in the taxonomy of Aspergillus section Nigri. Widely used β-tubulin primers are not specific for the benA gene for some taxa and preferentially amplify the tubC paralogue. Sequences of the tubC paralogue are widely combined with benA sequences in recent taxonomical works as well as other works, resulting in incongruent trees. In this study we newly provide benA sequences for several ex-type strains, which were characterised using the tubC gene only. We designed a highly specific forward primer to benA designated Ben2f for use in Aspergillus section Nigri, and tested specificity of numerous primer combinations to β-tubulin paralogs. The primer pairs with the highest specificity to the benA gene and functional across species in section Nigri includes Ben2f/Bt2b, Ben2f/T22 and T10/T22. We also provide tools based on codon usage bias analysis that reliably distinguish both paralogues. Exon/intron arrangement is the next distinctive characteristic, although this tool is not valid outside section Nigri. The species identity of taxa from the A. aculeatus clade used in previous molecular studies was revised using combined molecular data (ITS, benA, caM). These data together with two different PCR-fingerprinting methods indicated that A. japonicus should be treated as a synonym of A. violaceofuscus. Similarly, A. fijiensis is reduced to synonymy with A. brunneoviolaceus.


INTRODUCTION
Aspergillus section Nigri includes biotechnologically, economically and medically important organisms. Some species produce important organic acids or enzymes. Others are known to cause food spoilage or severe human infections. There are 26 species distinguished based on the polyphasic approach . A number of schemes have been proposed for classification and identification of black aspergilli. Molecular approaches have revealed high diversity among species that are difficult to recognise based solely on their phenotypic characters (Samson et al. 2007. The ITS sequence offers only low discrimination level to major clades. The benA and caM genes are the most informative loci in classification of black aspergilli and are broadly preferred by taxonomists and in routine identification. Despite this, the benA or caM genes alone are not able to distinguish all species in section Nigri . Thus, the polyphasic approach recommended for description of new Aspergillus species (Samson & Varga 2009) is only employed in a limited fashion in section Nigri.
Paralogous genes represent well-known problems in taxonomy that uses molecular features as the indispensable tool for delineation of taxa. The origins of paralogous genes are most frequently attributed to gene duplication. Among members of the tubulin superfamily, only α-, β-and γ-tubulins have homologues in fungal genomes (Dutcher 2001). β-tubulin paralogues are represented in Aspergillus by two genes designated benA and tubC. Gene benA (tub-2) encodes polypeptides designated as β1-and β2-tubulin (diversifying due to different post-transcriptional modifications) and is the third most utilised gene in fungal multilocus phylogenies (Feau et al. 2011). Polypeptide β3-tubulin (Weatherbee & Morris 1984) is encoded by the tubC paralogue and was first described in A. nidulans . Whereas benA is a housekeeping gene, tubC is not essential for growth and is expressed only under specific conditions (see Discussion). In this study, we show that β-tubulin paralogues are commonly mixed in taxonomical studies producing incongruent phylogenetic trees and bringing discrepancies into the taxonomy of black aspergilli.

Molecular studies
DNA was extracted from 7 d old colonies using the Microbial DNA Isolation Kit (Mo-Bio Laboratories, Inc.). For phylogenetic analysis, the ITS region, partial benA gene and partial caM gene were chosen because they were used in recent taxonomical monographs (Samson et al. 2004). The Mastercycler Gradient (Eppendorf) was used to amplify the desired regions. The ITS region of the rDNA was amplified using primers ITS1F and NL4 (O'Donnell 1992, Gardes & Bruns 1993). The partial caM gene sequences were amplified with primers CF1M or CF1L and CF4 as described by Peterson (2008). PCR amplification was performed using the ʿType Iʾ conditions (see below). The benA and tubC loci were amplified using primers listed below in section ʿPrimer specificity testingʾ and amplified using two different cycling conditions. The PCR product purification and sequencing was provided by Macrogen Europe, The Netherlands. DNA sequences obtained in this study were deposited in the EMBL database under the accession numbers listed in Table 1, while other sequences are deposited under numbers HE818079 -HE818087. β-tubulin paralogue tubC is frequently misidentified as the benA gene in Aspergillus section Nigri taxonomy: primer specificity testing and taxonomic consequences

Primer specificity testing
The Bt2a (5´-GGTAACCAAATCGGTGCTGCTTTC), T10 (5´-ACGATAGGTTCACCTCCAGAC) and Ben2f (5´-TCCA-GACTGGTCAGTGTGTAA) were selected as forward primers; Bt2b (5´-ACCCTCAGTGTAGTGACCCTTGGC), T224 (5´-GAGGGAACGACGGAGAAGGTG), T222 (5´-GACCGGG-GAAACGGAGACAGG) and T22 (5´-TCTGGATGTTGTTGG-GAATCC) were selected as reverse primers for testing. The position of all primers is indicated on Fig. 1. The specificity of all possible primer combinations was tested using two different PCR cycling conditions. The identity of amplification products was verified by sequencing or by length of fragments observed on the electrophoretograms. The mixture (25 μL) contained 50 ng of genomic DNA, 20 pmol of each primer, 0.2 mM of dNTPs, and 1 U of PerfectTaq DNA polymerase with the respective buffer. On an Eppendorf Mastercycler Gradient (Eppendorf, Hamburg), two types PCR cycling condition were used. The ʿType Iʾ conditions comprised 32 cycles under the following temperature regime: 95 °C/3 min, 55 °C/30 s, and 72 °C/1 min (1×), 95 °C/30 s, 55 °C/30 s, and 72 °C/1 min (30×) and 95 °C/30 s, 55 °C/30 s, and 72 °C/10 min (1×). In case of primer combinations with T22 primer, the extension time 90 s was used. The ʿType IIʾ touchdown cycling conditions involved an initial 2 min denaturation step at 93 °C, followed by 5 cycles in which the DNA samples were denatured at 93 °C for 30 s and annealed for 30 s with a decrease in 1 °C in each successive cycle. The regime started with an annealing temperature 65 °C decreasing to 60 °C. The extension was proceeded at 72 °C for 1 min (in case of primer combinations with T22 primer, the extension time 90 s was used). Annealing at the 60 °C was then used for further 33 cycles with a final extension for 10 min.
The phylogram is shown on Fig. 2 and includes all β-tubulin sequences of taxa belonging to the A. aculeatus clade and benA sequences of other species from section Nigri from GenBank. The introns were extracted from sequences using features for benA and tubC paralogues published by May et al. (1987) (reference annotated sequences: M17519 for benA and M17520 for tubC). The maximum likelihood (ML) method implemented in MEGA5 (Tamura et al. 2011) was used with settings similar to that for the phylograms in Fig. 3 (see below). In addition, the codon positions were indicated and all were included in the analysis. The tree with the highest log likelihood is shown (Fig. 2). There were a total of 207 positions (73 variables) in the final dataset.
The ITS region and partial benA and caM sequences were combined. There were a total of 1 346 positions in the final dataset (430 variables). The evolutionary history was inferred by using the ML. A discrete Gamma distribution was used to model evolutionary rate differences among sites. The product was identified as tubC based on the similar length of the product to the tubC from the ex-neotype culture of A. violaceofuscus (Fig. 3).
Sequences that originated from this study are in bold print. Fig. 1 The exon-intron arrangement of benA and tubC gene of A. aculeatus. The exons are in red (benA) and green (tubC) colour. The position of primers used for specificity testing (see Table 4) and their orientation is designated by arrows above and below exons.  NRRL 6362 was used as an outgroup. The tree with the highest log likelihood is shown (Fig. 3). The Bayesian tree inference analysis was used to calculate the posterior probabilities of branches (Huelsenbeck & Ronquist 2001). Separate partitions were created for exons and introns (benA and caM), ITS spacers and 5.8S rDNA. The exons were analysed allowing codon positions to be independent datasets. A general time reversible model was used with gamma-distributed rate variation across sites allowing six different types of substitutions. The MCMC analysis with 8 × 10 7 generations was run with two parallel chains incrementally heated by a temperature of 0.7, starting from a random tree. One tree was saved per 1 000 generations, and the run was ended when the likelihood scores of sampled trees approached similar values. A burn-in and convergence of the chains were determined with Tracer v1. Thermal cycling parameters were initial denaturing for 2 min at 96 °C, followed by 40 cycles of 30 s at 96 °C, 45 s at 44 °C, and 90 s at 72 °C, with a final elongation at 72 °C for 10 min. The amplified products were subjected to electrophoresis on 1.8 % agarose gels stained with ethidium bromide, and the banding patterns were visualised under ultraviolet light. The bands were scored as either present or absent for each strain. If more than two thirds of bands were shared between isolates, the patterns were evaluated as the same type and designated by the Roman numerals. If more than half and less than two thirds of bands was shared, the patterns were evaluated as subtypes and designated by lower case letter following Roman numerals of particular types.

Codon usage (CU) analysis
A dataset was used for CU analysis that involved all sequences used for construction of the phylogram in Fig. 2. The final alignment included 207 positions corresponding to 69 amino acid residues. CU characteristics were determined using program CodonW (http://codonw.sourceforge.net) (Peden 1999). The Codon Bias Index (CBI) (Bennetzen & Hall 1982) and frequency of optimal codon (F OP ) (Ikemura 1981) were calculated taking the optimal codons of A. nidulans (Lloyd & Sharp 1991) implemented in CodonW. Other characteristics calculated were the G+C composition in the third position of synonymous codons (GC 3s ) and the composition of particular bases in the third position of synonymous codons (A 3s , C 3s , G 3s , T 3s ).

Statistical analysis
Data analysis was performed with the software package PAST (Hammer et al. 2001). Data characterising CU did not show a normal distribution (based on Shapiro and Wilk's W-test), and therefore, the Mann-Whitney test (MWt) was used to test for significant differences in particular CU statistics between paralogues.

Source of isolates
The isolates used in this study were obtained from the Culture

TubC paralogues in the GenBank/EMBL/DDBJ database
Phylogenetic analysis (Fig. 2) and CU analysis (Fig. 4) indicated that β-tubulin sequences of black aspergilli cluster in two distinct clades corresponding to benA and tubC paralogues. Fragments delimited by the Bt2a/Bt2b primer pair corresponding to the tubC paralogue only have two introns (the intron following the 53rd amino acid residue from the beginning of both paralogues is missing in tubC) compared to three intron fragments in the benA gene (Fig. 1) as previously recognised by Peterson (2008). The affected species were A. aculeatus, A. brunneoviolaceus, A. fijiensis, A. japonicus and A. violaceofuscus. In addition, we amplified the tubC paralogue from one species that is tentatively new (Table 1, Fig. 3) using Bt2a/Bt2b primers. The species identities of some isolates for which only the tubC paralogue was published were questionable due to the lack of diagnostic caM and benA sequences. These taxa were sequenced and their revised identity is presented on Fig.  2, 3 and

F OP CBI
The sequences of the tubC paralogue are predominant sequences among β-tubulin sequences deposited in GenBank for A. japonicus and A. violaceofuscus (A. japonicus is here considered as synonymous to A. violaceofuscus -see below). For A. aculeatus, there was no benA sequence deposited that could be assigned to this species, although the species name ʽA. acu-leatusʼ was the most frequent name under which sequences of species from the A. aculeatus clade are deposited. Similarly for A. violaceofuscus, there were only two benA sequences belonging to non-type isolates; other sequences represented the tubC paralogue ( Table 2). The absence of benA sequences in databases for type specimens of A. aculeatus, A. japonicus and A. violaceofuscus is taxonomically important. Appropriate benA sequences were amplified in this study (see below) and deposited in the EMBL database ( Table 2).

Distribution of tubC paralogue in section Nigri
In section Nigri, the β-tubulin tubC paralogue is most likely only present in some taxa from the A. aculeatus clade sensu Varga et al. (2011). Concerning this clade, no tubC sequences were published for A. aculeatinus, A. indologenus and A. uvarum. PCR tests with A. aculeatinus F-596 and IHEM 20714, the only one of the three species mentioned and included in our study, confirmed the presence of tubC. The tubC paralogue was absent in A. carbonarius and A. niger, two taxa with complete genome sequences.

Taxonomical consequences
The sequences of the tubC paralogue used in combined datasets with benA sequences resulted in long, marginal, wellsupported branches in phylogenetic trees (Samson et al. 2004, de Vries et al. 2005, Varga et al. 2007, 2011, Noonim et al. 2008, Sørensen et al. 2011). The tree topologies were clearly different from those constructed based on caM sequences, and the relationships between taxa were distorted. Aside from taxonomical works, the tubC paralogue was also amplified in studies focused on medical or food mycology (Table 3).
Our data indicate that the position of A. japonicus and A. violaceofuscus as separate taxa that was proposed by Varga et al. (2011) is not supported by sequence data and by two fingerprinting methods used (Fig. 3). There are no unique positions in alignments shared across isolates that were designated as A. violaceofuscus by Varga et al. (2011) differentiating them from isolates of A. japonicus. A previous molecular study of Peterson (2008)  Based on the ITS and benA data, A. fijiensis is indistinguishable from A. aculeatinus and A. brunneoviolaceus (Fig. 3). The caM data separate A. fijiensis from A. aculeatinus and only one unique position in the caM locus segregates weakly supported clades with an ex-type isolate of A. brunneoviolaceus (CBS 621.78 T and IHEM 4062) and A. fijiensis (CBS 313.89 T , IHEM 22812 and NRRL 359) (Fig. 3). The intraspecies genetic distances between two A. fijiensis isolates (CBS 313.89 T and CBS 119.49) designated by Varga et al. (2011) are similar to those between ex-type isolates of A. fijiensis and A. brunneoviolaceus. Two fingerprinting methods previously used in Aspergillus (Nováková et al. 2012) and Penicillium (Tuthill 2004) for typification at species and subspecies level also showed no support for A. fijiensis as separate species (Fig. 3). Similarly, no unique morphological features differentiating A. fijiensis from A. brunneoviolaceus were found (Hubka, unpubl. data

BenA amplification in Aspergillus section Nigri and primer specificity
The primer combination of Bt2a and Bt2b is the most important in amplification of the benA gene in Aspergillus. Under standard conditions (annealing 55 °C), we found that this combination shows specificity for the tubC paralogue or both paralogues are amplified simultaneously (Fig. 5). After the annealing temperature is increased, the specificity for the tubC paralogue is increased (Fig. 6).
We tested all possible combinations of some previously published primers (Glass & Donaldson 1995, O'Donnell & Cigelnik 1997) marked on Fig. 1 and a newly designed primer Ben2f on a set of uniseriate as well as biseriate species from section Nigri (Table 4) using two different settings. The problems with primer specificity are completely solved when using the newly designed Ben2f primer as the forward primer (Fig. 5) in combination with Bt2b and T22 reverse primers (Table 4). Another primer combination functional across species from section Nigri and showing good benA specificity was the T10/T22 combination.

Codon usage
There were significant differences in codon usage bias between paralogues benA and tubC. CBI and F OP statistics were significantly higher (MWt p-values < 10 -13 ) for the benA gene (Fig. 4) indicating its higher level of codon bias. Other statistics that were significantly different between paralogues comprised GC 3s , A 3s , G 3s , T 3s (all MWt p-values < 10 -13 ).

Paralogous genes and taxonomy
The combination of sequences belonging to paralogous genes with non-homologous functions in the same phylogenetic analysis is a great risk and might cause incongruences within and between datasets. Taxonomists prevent the impact of paralogous genes on taxonomic conclusions by using a polyphasic approach. Molecular data from several non-linked loci are combined with morphological, physiological and other traits to define interspecies boundaries (Samson & Varga 2009). Regarding quickly speciating species such as Aspergillus (sequence divergences between sibling species are mostly 0-5 % depending on the locus examined), the importance of morphology is often only secondary. Together with the growing number of cryptic species, taxonomy misses an important tool for elimination of the impact of paralogous genes on taxonomical conclusions.
In fungal genomes, an inconstant number of β-tubulin paralogues can be found (Hubka 2011) that can be randomly amplified when using primers with low specificity. Keeling et al. (2000) tried to construct a phylogeny of Fungi based on β-tubulin gene sequences. In several cases, the authors amplified two or three paralogues for some species using newly designed primers. The results of such an analysis and its interpretation had a very limited value. β-tubulin primers for Fungi were also designed by Einax & Voigt (2003), but their specificity is disputable due to the number of by-products on depicted electrophoretograms. Peterson (2008) first noted differences in intron numbers between β-tubulin amplicons of uniseriate black aspergilli amplified by the Bt2a/Bt2b primer pair. Because of doubt about the homology, the sequences of β-tubulin genes were not used in multilocus analyses of sections Nigri, Usti and Nidulantes. This finding was not further kept in mind by taxonomists, and fragments with a variable intron number were combined in recent taxonomical studies (Table 3). Peterson´s finding was misinterpreted by Sørensen et al. (2011) as a variation in intron number in the benA gene of A. aculeatus and A. violaceofuscus, compared to the other members of section Nigri. In fact, there is no difference in the number of introns between black aspergilli for the benA gene. Most of the Aspergillus species produce three intron Bt2a/Bt2b benA fragments. The homologous fragments of the tubC paralogue in section Nigri include only two introns (Fig. 1). The differences in fragment lengths can also be observed on electrophoretograms (Fig. 5, 6). Nevertheless, intron number cannot be used as reliable marker for distinguishing of both paralogues as presumed by Peterson (2008). Across Aspergillus species fragments can be found from benA as well as tubC with two or three introns in fragments bordered by the Bt2a and Bt2b primers (Hubka & Kolařík, unpubl. data I   II   I   II   I   II   I   II   I   II   I 3   II   I   II   I   II   II   I   II  Two different PCR cycling conditions were used for amplification; Type I -ʿstandardʾ PCR conditions, Type II -touchdown cycling conditions; for details see Material and methods.
3 Amplification products are relatively weak, it can be improved by increasing the number of PCR cycles to 40, although tubC fragment can occur in A. violaceofuscus and A. aculeatus.
Amplification products: b = benA; t = tubC; bt = benA and tubC fragment together; mp = multiple non-specific products; 0 = no product. Table 4 Selected primer combination and their specifity to β-tubulin paralogues across species belonging to Aspergillus section Nigri.
specificity in some taxa. Decreased divergence between paralogues in taxonomically important regions bordered by Bt2a/ Bt2b primers in section Nigri most likely participated in altered specificity of the Bt2a/Bt2b pair that preferentially amplifies tubC paralogues or both paralogues at the same time. We solved this problem with non-specifics primers in section Nigri by implementation of a new Ben2f primer that showed good benA specificity and wide usability across black aspergilli (Fig. 5, Table 4). Further studies are needed to verify the specificity of β-tubulin primers in other sections of the genus Aspergillus and also in other fungi. We can also hypothesise that taxonomy of other fungal groups may be affected by illegitimate use of paralogous genes analogically to black aspergilli. The GenBank database includes more than 25 000 β-tubulin sequences, and we only examined a subtle fraction.
High divergence of β-tubulin paralogues was first observed in A. nidulans by May et al. (1987). This divergence is accompanied by a different codon spectrum used by paralogues. The benA gene is highly biased and its codon spectrum is markedly limited. In contrast, no preferences in the use of synonymous codons can be observed in tubC. As we show here, this different level of codon bias can be used as an excellent marker for distinguishing benA and tubC. CBI and F OP statistics characterise the level of codon usage bias by one value and are suitable for gene comparison (Fig. 4). In addition, optimisation of both characteristics is available for A. nidulans (Peden 1999). F OP is a simple ratio between the frequency of optimal codons (that appear to be translationally optimal) and the total number of synonymous codons. It ranges from 0 to 1 (when a gene is entirely composed of optimal codons). CBI is a measure of codon bias towards a subset of optimal codons and is similar to F OP . In a gene with extreme codon bias, CBI may equal 1.
Codon usage bias also has important functional consequences.
Highly biased genes are generally highly expressed and perform important functions in contrast to genes with low levels of codon bias (Sharp et al. 1986, Sharp & Devine 1989. It was also demonstrated that genes with similar codon usage are usually co-expressed during the life cycle (Lavner & Kotlar 2005, Najafabadi et al. 2009). This is in agreement with previous observations regarding functions of β-tubulin species. Products of the benA gene play an important role during whole vegetative growth: they participate in the formation of the mitotic spindle  Table 1 and designated A-L. The specificity of the Bt2a and Bt2b primer pair is apparently low in contrast to primer pair Ben2f and Bt2b that is highly benA specific.  and in movement of organelles including the nucleus (Oakley & Morris 1980, 1981. The product of the tubC paralogue participates in conidiogenesis in A. nidulans but is not essential for this process , May 1989. Further studies are needed to uncover if notably different divergence between paralogues in A. nidulans and in black aspergilli also has functional consequences. Searching for sequence similarity via Blast (http://blast.ncbi. nlm.nih.gov/Blast.cgi) offers a very simple method for discrimination between both paralogues. Annotated sequences labelled as benA or tubC that were amplified in this study are available under the accession numbers listed in Table 1 and 2. The divergence between benA and tubC is sufficiently high for their clear differentiation.

Taxonomic notes
Aspergillus aculeatus and ʿA. japonicusʾ are well-supported species based on molecular data (Pařenicová et al. 2001, Samson et al. 2007, although both species are indistinguishable based on morphology (Hamari et al. 1997). In the past, A. violaceofuscus was treated as a colour variant and synonym of A. aculeatus by Raper & Fennell (1965). Varga et al. (2011) treated A. violaceofuscus as a valid species related to A. japonicus based on a polyphasic approach, although the analysis of the β-tubulin locus was based on a mixed benA/tubC dataset. We re-examined the benA, tubC and caM sequences of isolates treated as A. violaceofuscus and A. japonicus by , and we did not observe any molecular support for separation of A. japonicus and A. violaceofuscus. Although differences in conidial shape can be observed in the neotype culture of A. violaceofuscus CBS 123.27 and the ex-type culture of A. japonicus (CBS 114.51), unique DNA characters at multiple loci should be present as a gold standard for Aspergillus species delimitation (Samson & Varga 2009). Additionally, the ex-type cultures share the same banding patterns provided by fingerprinting methods (Fig. 3).
Aspergillus brunneoviolaceus (CBS 621.78 = NRRL 4912), described by Batista & Maia (1955), was treated as a synonym of A. japonicus by Raper & Fennell (1965). We examined sequence data provided by Peterson (2008) that indicated that A. brunneoviolaceus should be a valid species. This taxon was omitted by Varga et al. (2011), and two very closely related isolates were proposed as a new species; A. fijiensis, although there is very low phylogenetic support for its delimitation (Fig. 3) from A. brunneoviolaceus. The extrolite data do not clearly support A. fijiensis as a separated species. Secondary metabolites produced by A. fijiensis differ among isolates CBS 119.49 and 313.89 that were used for description (Pařenicová et al. 2001). Additional differences are found in an ex-type isolate of A. brunneoviolaceus (Pařenicová et al. 2001). Consequently, A. fijiensis is here treated as synonymous to A. brun neoviolaceus. The position of another related isolate CBS 620.78 remains unresolved (Fig. 3).
Although these taxonomic remarks are not in direct consequence with illegitimate use of the tubC paralogue, confusion associated with tubC was introduced in Aspergillus section Nigri taxonomy and complicates description of new uniseriate taxa. There is a substantial call for searching for new molecular and physiological markers that are usable in the classification of black aspergilli. The concept of several recently described species is only based on genetic differences at one locus. The morphological concept is insufficient and the extrolite data are not fully resolved due to a number of newly described or revived species in the A. aculeatus clade, and intraspecific differences are found depending on the isolate tested and the methodology used (