In the centres of diversity, on the other hand, crop evolution is an ongoing phenomenon: the "loss" and origin of "new" alleles frequently occur; indeed, land race cultivar dynamics are quite high [149,150]. Even in extensive collections from the Mediterranean region, researchers have found low geographic coverage and poor representation of on-farm genetic diversity in ex situ collections established in Europe [151]. This is perhaps because as for wild species diversity, a stratified sampling strategy is required for full coverage [152]. Notwithstanding, no matter how well sampled, ex situ collections do not represent the natural population in the diversity hotspot, after a few decades. Therefore, a regular and systematic monitoring programme of the landrace cultivar population dynamics in well-defined in situ hotspots is needed to better understand the drivers of change [150] and to supplement existing collections. Systematic monitoring programmes would enable proper sampling at regular intervals to capture the changing allelic frequencies, which is particularly important for meeting breeding challenges in a changing climate.

#### 6.3.2. Identification of Duplicates

With large numbers of accessions in genebanks, another issue confronted by curators is the duplication of samples resulting from poor passport data and collection strategies, inconsistent documentation, and a lack of characterisation, amongst other factors. Labelling errors, hybridisation between stocks, and confusion about origins have also been identified as problems in ex situ collections. With over 7 million existing accessions in genebanks around the world [145] and increasing storage demands and costs, methods for efficient characterization and curation are required to avoid duplication. Genebanks have used molecular markers, such as simple sequence repeat (SSR), single-nucleotide polymorphism (SNP), amplified fragment length polymorphic (AFLP), and chloroplast DNA markers to characterise collections for their diversity and to identify duplicates with the aim of rationalising future collection and conservation efforts.

Lately, genome sequencing is being used to identify new patterns of variation [153,154], alleles of interest for breeding programmes [155], as well as duplicates [156] in collections. Methods are being developed to sequence large populations at low cost, including complex polyploid genomes [157–159]. Additionally, genome sequencing enables the unravelling of patterns of evolution of ancient crops, such as apple [160], grape [161], rice [114], and wheat [162], giving insights to their progenitors and conservation needs. Furthermore, they enable quantitative trait mapping and novel allele mining from large genetic collections [163]. When applied across different genebanks, these methods will enable the cost-effective and efficient management of germplasm and better stewardship of valuable genetic resources.

#### 6.3.3. What to Conserve and Use—The Concept of Core Collections

It is practically impossible to cover the entire range of landraces or the diversity within CWR in a genebank. To manage the increasing number of accessions in collections and the resultant difficulties in monitoring, regeneration, evaluation, etc., the concept of core collections arose [164]. Identifying 'representative' samples within large collections helps to better utilise genetic resources when there are tens of thousands of accessions to choose from with just basic passport data. Such subsets in a collection represent the maximum diversity without redundancy. This allows the supply of a set of accessions for evaluation or breeding purposes without compromising the diversity within the collection.

With the introduction of the core set concept, genebanks now has a tool to manage their collections more efficiently. Core collections have since become focal points for conservation prioritisation, phenotypic evaluation, genotyping, and exchange. In the early years, the objective was to develop a single entry point for users, providing them with the widest diversity within a manageable number of accessions. Collections were defined hierarchically using taxonomic characterisation, genomic distance, and geographic data, dividing the collection into clusters and selecting 'core' genotypes within those clusters using methods such as proportional allocation, log frequency allocation, and constant allocation groups. Once a core collection is developed, it can be validated using several methods, such as mean comparisons with the entire collection, the homogeneity of variances and frequency distributions among traits between the entire collection and the core, and optimising correlations [165].

Alternatively, diversity can be assessed using molecular markers, and core sets developed using a maximization of the number of alleles observed in each marker locus without relying on a stratification strategy. This M (maximization) strategy examines all the accessions for the alleles in the collection and identifies those that maximize the number of observed alleles at the marker loci. These can then be chosen as final candidates for the core. The superiority of this marker-based method is derived from the correlation between observed allelic richness at the marker loci and allelic richness at other loci. The software uses iterative procedures to select samples with the highest diversity as measured by the number of alleles and the trait classes that account for the greatest proportion of the collection variability based on the M strategy. The method can also be applied for both quantitative and qualitative data [166]. For example, this method enabled the production of a mini core subset of 217 accessions representing a core of 1794 accessions from the United States Department of Agriculture Agricultural Research Service (USDA-ARS) collection of

rice, comprising more than 18,000 accessions [167], and to produce a core of 20 accessions within 450 apple accessions of the New Zealand apple germplasm repository, which was targeted for the first round of cryopreservation using winter dormant buds [168,169].

Core collections can also be developed using a multivariate distance approach [170] to develop the clusters followed by selection within those. Applying these methods, the International Crops Research Institute for the Semi-Arid Tropics (ICRISAT) has developed core sets from global collections of eleven crop species in ICRISAT genebank accessions [171]. Subsequently, mini core collections have been developed within the core collection using the same principles.

Core subsets for global collections, such as those maintained by the CGIAR, often serve as reference panels so that researchers around the world can evaluate the same genetic resources in different environments. Individuals in those reference sets may be selected as "controls" for phenotyping efforts (particularly relevant for disease and pathogen resistance/susceptibility and adaptation to different agro-ecological zones) so that results can be compared among research and breeding programmes. Some plant collections, particularly those that are clonally propagated, are difficult to distribute across international borders. Hence, multinational plantings of international core subsets for clonally propagated collections would help ensure that international communities have access to the same plant materials for comparative research.

#### *6.4. Types of Collections in Genebanks, Their Management and Utilization*

Ex situ collections are conserved in different forms depending on the reproductive behaviour, genetic composition of the population, use of the crop (seed, fruit, fibre, vegetative materials, etc.) and other considerations. These aspects are briefly discussed below, with references to recent studies for details.

#### 6.4.1. DNA Banks

Plant DNA banks were initially developed to create genetic libraries for evolutionary studies, to understand biological diversity, and to collect genomic information [172]. However, with habitat loss, species extinction is happening at a rapid rate, particularly in areas with high genetic diversity. Therefore, DNA banks are increasingly used to store genomic and diversity information of species in these vulnerable hot spots. Additionally, our ability to extract DNA from fossilised plant remains is providing new research opportunities in paleoecology, phylogeography, and evolution, including crop domestication. Therefore, the number of DNA specimens from extinct species is increasing in DNA banks.

Genomic research has seen unprecedented advances in recent years. However, the physical DNA from the published sequences is generally not accessible to researchers. Access to the original samples is important for conducting new studies, to extend or complement existing results, and to support good scientific practice, enabling the verification of published results. To address these needs, the Global Genome Biodiversity Network (GGBN) was formed in 2011, with the aim of developing high-quality, well-documented, and vouchered collections that store DNA or tissue samples, as well as to encourage and enable scientists to complete documentation chains between vouchers, tissues, physical DNA, sequences, and publications [173,174]. The GGBN currently has 99 member organisations from 35 countries, with over 4.1 million samples from 5152 families, covering 38,074 genera and 140,182 species. Their updates can be viewed in annual newsletters [175] and the data can be accessed through their data portal [176].

DNA banking is considered the most economical way of conserving genetic information of plant genetic resources, and is also the easiest way of exchanging genetic information across borders because DNA samples occupy less space, are more stable, and phytosanitary certification is not needed.

#### 6.4.2. Orthodox Seeds

There are three types of seeds from a conservation perspective. Orthodox seeds can be dried to low moisture contents (about 3% for oily seeds and about 7% for starchy seeds) without damage and be stored dry at low temperatures without losing their viability over long periods. Recalcitrant seeds cannot withstand desiccation to moisture contents below 20%. Intermediate seeds can be dried to a moisture content of 10% to 12%, but further desiccation reduces viability and/or dry seeds are injured by low temperatures. Seeds of most of the main crop species belong to the orthodox category, meaning that they can be safely cooled to standard long-term storage conditions of −18 ◦C without losing viability once dried to appropriate moisture contents. Seed storage is well researched and is the most efficient and cost-effective method of the conservation of plant germplasm, with about 96% of accessions held in genebanks worldwide as seeds [177]. In most of the genebanks, the seed vaults are maintained at −18 ◦C for long-term storage, whereas active collections are maintained at 5 ◦C to 10 ◦C. Seeds may be stored in dedicated cold rooms (Figure 8) or in domestic deep chest freezers or refrigerators, for which genebank standards and technical guidelines have been published by the FAO [178].

**Figure 8.** Seed vaults for long-term storage at the United States Department of Agriculture—Agricultural Research Service genebank in Fort Collins, Colorado.

There are many advantages of conserving germplasm as seeds, including better security than in the natural environment, less space required, methods being simple, easily accessed and exchanged, and long storage periods. Additionally, seed collections capture more allelic diversity of the population than clonal collections. Therefore, in the case of fruit crop wild relatives, seed conservation is an option. However, as mentioned above, the evolutionary processes are frozen in time: as in any ex situ conservation method, regeneration is required from time to time. There can be gaps in the collection and initial set-up, and ongoing maintenance can be costly, including providing a constant power supply.

Under the ITPGRFA multilateral system, there are over 730,000 accessions available from the collections of the CGIAR system. The majority of these accessions are held in the form of seeds, with only 23,862 conserved as clones in vitro and 29,122 in field collections [143]. The number of accessions according to the species and the 11 participating genebanks are listed [143]. Additionally, the IBPGR has also coordinated the creation of an internationally linked system of genebanks known as the Registry of Base Collections (RBC) to conserve and make a subset of those materials available. Under this system, 144,000 accessions are available in 52 selected genebanks spread across all continents, covering 80 genera and 250 species [179]. Engels and Ebert [177] recently critically reviewed the current status of global seed banking with recommendations for improvement, emphasising the role of functional genomics and phenomics as well as strengths and weaknesses within regulatory frameworks and the strategies for linking national programmes with the global network.

#### 6.4.3. Genebanking of Clonal Crops and Recalcitrant Species Field Collections

The traditional approach for maintaining germplasm of clonal species, species that are sterile or semi-sterile, and those producing recalcitrant seeds, is in field genebanks. They provide an excellent opportunity for curators to assess the diversity by phenotyping and for plant breeders to directly use accessions in their breeding programmes. Woody perennial crops (WPCs) represent about one-half of crop plant diversity and one-third of the 167 major crops. WPCs are usually clonally propagated because they are obligately outcrossing. Therefore, the conservation of germplasm through seed banks is not applicable to maintain the genetic characteristics of different heterozygous individuals. Moreover, the seeds of many species are recalcitrant, and therefore cannot be stored in traditional seed banks. Additionally, WPC species have a juvenile stage that can last for several years. In fact, plants generated from seeds usually exhibit strong, undesirable juvenile characters, such as a thorny habit and delayed flowering and fruit production [180]. In order to avoid these problems, the long-term conservation of WPC germplasm can be effectuated as field collections with clonal material in order to maintain the elite genotypes that form the foundation of woody perennial agriculture. However, WPCs represents only 5.8% of ex situ germplasm collections. Despite their importance, field collections are expensive to establish, and maintenance requires high costs for specialized technical personnel and land. Usually, woody perennials have a large plant size and therefore need large areas for maintenance in the field. Moreover, field collections have high risks of loss because they are exposed to natural disasters and are subject to biotic and abiotic stresses. In particular, the risk of pathogen infections transmitted through vegetative propagation is high and difficult to avoid. There are many examples of such cases, e.g., the loss of accessions in the apple germplasm collection in USDA due to fireblight (*Erwinia amylovora*) [181] and the loss of kiwifruit germplasm in New Zealand due to the incursion of *Psuedomonas syringae* pv. *actinidiae* [182]. Therefore, many genebanks back up their field collections in separate locations and also, more reliably, using lab-based conservation methodologies such as in vitro storage or cryostorage.

Establishing field genebanks may seem straightforward, but there are established best practices for sampling, the duplication of collections, and cataloguing the accessions. These can be found, for example, in the IPGRI training manual on "Establishment and management of field genebanks" [183], which has separate chapters on principles, legal issues, plant health, choice of materials, genetic considerations, planting layout, management characterisation and evaluation, utilisation, and economics. Another useful document is the IPGRI handbook No 7 on "Technical guidelines for the management of field and in vitro genebanks" [184]. Field genebanks of horticultural crops are a long-term commitment. Therefore, in many countries, field collections and genebank operations are undertaken by government departments.

#### In Vitro Collections

In general, lab-based germplasm conservation strategies fall under two major categories: slow growth procedures and cryopreservation. Both strategies require efficient regeneration system with high regeneration efficiency via organogenesis and/or somatic embryogenesis. In vitro storage has the advantage of maintaining collections under diseasefree conditions and the cultures have fewer biosecurity requirements when material is exchanged across borders. On the other hand, because collections have high genetic diversity, the response to tissue culture media will be variable, which therefore requires prior research to optimise media for different species and even genotypes within species. There have been many studies on the tissue culture of cultivated species; therefore, this is not as daunting as it used to be a few decades ago. In vitro genebanks enable the rapid multiplication of material when required and provide a safe environment for managing germplasm collections in a confined space away from the field. To avoid somaclonal variation (genetic variation induced under tissue culture conditions), pathways using dedifferentiation and

adventitious regeneration should be avoided. Therefore, intact shoot tips and axillary buds are typically used with minimal use of plant growth regulators.

After initiation and the multiplication of accessions in tissue culture, they are stored under slow-growth conditions to increase the period between subcultures, thus significantly reducing the cost of labour and materials. Generally, this is achieved by a combination of several factors: (a) decreased light (generally about 5–10% of the standard culture conditions); (b) low temperature between 4 and 21 ◦C with tropical species requiring higher temperatures; (c) lower concentration of mineral nutrients and sucrose in media, often without growth regulators; (d) smaller size of culture vessel/vial; (e) inclusion of growth retardants in media; and (f) osmotic stress using chemicals [62,185,186]. Unless available in the literature, these conditions need to be determined for each species by experimentation. Many genebanks around the world have in vitro storage for some species. Genebanks holding large numbers of accessions in vitro are given in Table 4.

Synthetic seed technology (SST) is an innovative and sustainable approach to preserve the biodiversity of clonally propagated woody perennials. A synthetic seed is defined as an artificially encapsulated somatic embryo, vegetative bud, or any other micropropagule that plays the role of a seed and has the ability to give rise to a complete true-to-type plant. SST offers opportunities to conserve clonal genetic resources safely at low cost.

General guidelines for the storage [187] and status of in vitro storage under the CGIAR network, particularly CIP, IITA, CIAT, and Bioversity International [188], have been published. Additionally, in vitro techniques have proven useful in collecting germplasm when seeds are not available (off-season), or in situations where seeds are not likely to remain viable because of their recalcitrant nature. Pence et al. [189] recently reviewed basic approaches and principles of in vitro collection, and one of the crop species for which in vitro collection is routinely used is coconut.

**Table 4.** Genebanks with large in vitro collections.

Cryopreserved Collections—Stopping the Biological Clock

Cryopreservation is the storage of biological samples in liquid nitrogen (LN) at −196 ◦C or in its vapour phase (LNV) at −165 ◦C to −196 ◦C. The demonstration that winter shoots of woody species can be conserved in liquid nitrogen in 1960 by Sakai [197] gave the impetus to undertake research on cryopreservation of fruit tree germplasm [198]. In the meantime, Latta [199] demonstrated the survival of carrot and sweet potato cell cultures in LN when pre-treated with high concentrations of sucrose followed by freezing to −40 ◦C and subsequent transfer to LN. This method came to be called classical two-step freezing or slow freezing. However, it was the advent of vitrification methods in the 1990s that allowed genebanks to take up cryopreservation on a mass scale as it became applicable across different families and genera. The method relies on the dehydration of cellular content to the extent that sudden freezing does not allow water molecules to form lethal ice crystals, but enters a state of metastable glass—hence the term vitrification. Vitrification can be applied to naked meristematic explants, such as shoot tips and embryogenic cells, or to propagules protected by encapsulation in alginate beads [200,201]. The most widely used vitrification solution is Plant Vitrification Solution 2 [202] and the method is droplet vitrification, where propagules are held on an aluminium strip covered in a droplet of vitrification solution and directly immersed in LN. However, the recently developed V cryo-plate technique, where droplet vitrification and encapsulation techniques are combined, has shown an improvement in recovery over droplet vitrification in some crops [203–205]. The two methods as applied to grapevine have been described by Bettoni et al. [206]. Another advantage of cryopreservation using vitrification is its ability to eliminate viruses, phytoplasma, and bacteria [169,207–213] for the delivery of high-health propagation material.

Through cryopreservation, viable explants can be brought to a state where cellular division and metabolic processes are minimized to the extent that they cease, preserving the structure and function of the biological system—virtually stopping the biological clock. There are no biochemical processes or gene expression; therefore, the genetic material is safe from any changes, and is hence ideal for conservation. Most of the protocols depend on a tissue culture phase (except the dormant bud technique used mainly for apple cryopreservation); thus, it is important to ensure the precautions mentioned in the previous section are taken to avoid somaclonal variation. It is also cost-effective to maintain collections for extended periods of time in LN compared with field or in vitro collections [214], and the cost effectiveness increases as more accessions are added to the collection [215]. LN is used to freeze material, so the method is not dependent on an electric power supply; hence, it is an attractive method for countries with an insecure power supply, although it is dependent on a reliable LN supply.

The use of cryobanking in genebanks using vitrification methods has recently been summarised by Wang et al. [200], and for the USDA, currently the world's largest cryo collection, by Jenderek and Reed [216]. Panis [217] recently summarised the major cryopreserved collections worldwide. It is estimated that about 10,000 accessions are in long-term cryostorage using explants from in vitro culture material (mainly vitrification methods). Of these, over 80% belong to five species: potato (38%), cassava (22%), bananas and plantains (11%), mulberry (12%,) and garlic (5%). Other large collections are in apple, using winter dormant buds. Cryopreservation is also used for intermediate and recalcitrant seed crop conservation, such as coconut [218] and coffee [215], as well as for pollen. Pollen is naturally dehydration-tolerant and can be used to conserve the nuclear genetic diversity of CWR, recalcitrant species, endangered and rare species, and fruit and ornamental crops. The cryopreservation of pollen allows access to pollen when needed by breeders, particularly useful when breeding lines have differing flowering times, the use of CWR in breeding programmes, and for hybrid seed production programmes. Cryopreserving the pollen of male sterile lines for use on female lines can save large amounts of land in hybrid seed production fields and labour to collect pollen during busy periods of hybridisation. Pollen of many tropical plant species has successfully been cryopreserved: for example, the Indian Institute of Horticultural Research in Bangalore holds 650 samples of pollen from 45 species belonging to 15 plant families [219].

Synthetic seeds can be stored for long periods using vitrification-based cryopreservation [220,221]; both somatic embryos and other somatic tissues with meristematic regions

have been used in cryopreservation and methods have been optimized to achieve post-thaw regeneration rates that meet genebank standards for the implementation of cryopreservation [201,206,222,223]. However, high genetic variability in the somatic embryogenesis response [224,225] is a barrier to use this propagule in cryopreservation.

Cryopreservation is only applied to a limited number of crops in some tropical genebanks outside the CGIAR system. Interestingly, the banana genebank, including the cryo-collection, is in Europe [226]. Since the introduction of vitrification methods, cryopreservation has become operationally simple and easily adaptable to any laboratory with tissue culture facilities. Collections of many tropical species, such as cassava, banana, sweet potato, and taro, are already being cryopreserved for their long-term security. The challenge is with CWR plus rare and endangered species of the tropics because their seeds are often recalcitrant and therefore cannot be preserved in traditional seed banks. Botanic gardens may only have one or a few specimens of each species in their field genebanks; however, the cryopreservation of seeds/embryonic axes and pollen enables the conservation of a much broader genetic diversity of CWRs or endangered species. Genebanks therefore need to invest in infrastructure and human resources for cryo-conservation. The integrated conservation strategies described above were recently used to save iconic New Zealand Myrtaceae species after the incursion of myrtle rust (*Austropuccinia psidii*) into New Zealand in 2017 [227]. The techniques used included the cryopreservation of pollen of *Metrosideros excelsa* and zygotic embryos of recalcitrant *Syzygium maire*, along with the successful hand pollination of flowers of *Metrosideros bartlettii*, a critically endangered species with only a few plants in the wild. Ex situ conservation strategies in field genebanks, in vitro slow growth, and the advantages and disadvantages of cryopreservation have recently been reviewed by Panis et al. [226].

#### **7. Challenges to Plant Breeding in Search of Right Germplasm**

Germplasm resources including CWRs carry alleles that are integral to creating new crop cultivars that can meet increasing consumer and environment demands. As changes in environmental conditions accelerate, so does the need for germplasm resources to breed more environmentally resilient crops. Properly characterised germplasm with information available to the plant breeding community will enable the deployment of improved cultivars on an ongoing basis. New cultivars with better water use efficiency, tolerance to soil toxicities, fertiliser response, and pest and disease tolerance will enable the replacement of at least some part of the costs for fertiliser, irrigation, and pesticides. Selection of the right genotype is the most effective means to achieve yield and quality improvements in crops, keeping other inputs to a minimum. Although plant breeders are often well aware of the necessity of maintaining genetic diversity in their breeding populations, they may lack the information to determine which of the thousands of accessions of a given crop would prove most beneficial for their breeding objectives. They may also be reluctant to introduce unadapted germplasm, with potentially negative impacts on quality, into their elite breeding materials. In some cases, they may lack the technical expertise or facilities to make interspecific crosses, for example, between CWR and cultivars of different ploidy levels that often require embryo rescue, in vitro pollination, and other interventions. Therefore, plant breeders should have access to facilities and expertise in plant biotechnology for the better utilisation of PGR.

Many international genebanks have a mandate to supply germplasm if the requests are justified. With more and more accessions being characterised and accessible databases made available, plant breeders should be trained in accessing the germplasm that matches their needs. The USDA-ARS holds over 576,000 accessions from 15,116 species and in 2015 alone, distributed over 239,000 accessions on request to national and international researchers [228]. The accessions available for supply can be searched in the Germplasm Resources Information Network-Global (GRIN-Global) database and ordered [229]. The CGIAR system of 11 genebanks holds over 736,000 accessions, with wheat and rice accounting for >196,000 and >144,000, respectively. Other major collections are in sorghum

(>39,200), beans (~38,000), barley (>31,500), maize (>29,700), pearl millet (>23,000), and forages and fodder (>44,000). The CGIAR system distributed over 3.9 million samples over a 10-year period from 2007 to 2016 [143]. The CGIAR system of networks has developed a Global Information System (GLIS), with an emphasis of assigning a unique digital object identifier (DOI) to each accession; these are linked to the GRIN-Global searchable database and to Genesys—an online database for global plant collections managed by the Crop Trust [230]. The VIR collection holds 346,666 accessions of PGR and CWRs, representing 64 botanical families, 376 genera, and 2169 species. Annually, 12,000–14,000 accessions are distributed by the genebank, of which 2000–3000 accessions are supplied outside Russia [25].

#### **8. Conclusions**

Significant improvements have been made for the better management of PGRs within international and national genebanks in the last few decades. More attention is being paid to managing 'difficult' species with recalcitrant seeds, clonal crops, and CWRs. Effective research and development mechanisms and policies have been established for the protection and conservation of biodiversity worldwide. Many countries have signed international treaties on biodiversity, as well as PGR conservation and exchange. Coordination among genebank curators, breeders, farming communities, government organisations, and different stakeholders should be strengthened to meet global obligations for sustainable management and the use of PGR for food and nutritional security. Genebank accessions are of no use unless they are accessible online or can be observed in the field. Therefore, conserved PGRs need to be characterised using advanced technologies and results with passport data made available for online access by users. The emphasis should be on duplicating collections for better security, creating and securing core collections, identifying gaps in collections and remediating them, and using advanced and integrated strategies for the conservation and dissemination of information. CWR and landrace varieties need to be secured for the future, using multiple strategies including in situ, on-farm, and ex situ conservation. Further exploration and collection of PGR from diversity-rich centres should also assume priority, recognising the high rate of loss of genetic diversity.

**Author Contributions:** Conceptualization, R.P. and F.C.; writing-original draft preparation, R.P. and F.C.; Writing-review and editing, R.P. and F.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

### *Article* **Seasonal Xylem Sap Acidification Is Governed by Tree Phenology, Temperature and Elevation of Growing Site**

**Manuel Pramsohler 1,\*, Edith Lichtenberger <sup>2</sup> and Gilbert Neuner 2,\***

**Abstract:** pH of xylem sap (pHx) was determined in three trees (*Malus domestica* (apple tree)*, Picea abies* and *Pinus cembra*) in response to seasonal changes. Conifer trees from lowland (600 m) were compared to trees growing at the alpine timberline (1950 m a.s.l.). Xylem sap was extracted with a Scholander pressure bomb and pHx was measured with a pH microsensor. In all species, pHx changed markedly with season. In spring, pHx was acidic; during winter, the pHx was more alkaline. In apple trees, the pHx did not show a significant correlation with temperature but was rather affected by developmental stage. During flushing in spring, xylem sap acidification took place concomitant to the developmental stage "tight cluster", when foliar development enables a significant transpiration and a consequent movement of water in the xylem. The xylem sap of the two studied conifers showed a significantly larger seasonal alkalinisation (+2.1) than found in apple trees (+1.2) and was significantly more pronounced at the timberline. Xylem sap acidification took place before bud break. pHx had a significant negative correlation with soil temperatures and corresponded to already reported pHx of angiosperms. Overall, pHx appears to be a sensitive stress marker and indicator of activity status in tree xylem.

**Keywords:** acidification; alkalinisation; bud burst; freezing; *Malus domestica*; pH; *Picea abies*; *Pinus cembra*

#### **1. Introduction**

The pH of the xylem sap (pHx) of plants ranges from acidic values of 4.5 to about 7.4 [1], which is in contrast to the milieu of the symplast, which has an alkaline pH ranging from 7.2 to 7.5 [2]. Only a few studies report on seasonal changes of pHx. In deciduous temperate trees, xylem sap usually becomes more alkaline in winter and then acidifies during spring [3–9]. In contrast, it can also be the other way round, as in *Juglans regia* pHx acidifies during winter and becomes more alkaline in spring [10].

Various mechanisms are reported to be responsible for the regulation of the pHx. All factors that affect proton pump activity [8,11,12] and xylem sap composition, especially the concentration of cations/anions and organic compounds, might be involved [13–15]. In addition, the pHx influences the concentration of dissolved CO2 in the xylem sap and therefore regulates the amount of CO2 that can be transported through the xylem in trees [1,9,16]. Seasonal pHx changes in deciduous species may be linked to the onset of the ascent of xylem sap, which is coupled to transpiration in the foliated stage. Unfortunately, studies reporting on seasonal changes of pHx in deciduous trees lack information on phenology and short-term variations during flushing in spring. For evergreen trees and gymnosperms, where the end of winter dormancy and the ascent of xylem sap are not linked to the formation of new leaves, only very few reports about seasonal changes of pHx are currently available. In conifers growing at the alpine timberline a marked pHx alkalinisation during winter is reported [17].

**Citation:** Pramsohler, M.; Lichtenberger, E.; Neuner, G. Seasonal Xylem Sap Acidification Is Governed by Tree Phenology, Temperature and Elevation of Growing Site. *Plants* **2022**, *11*, 2058. https://doi.org/10.3390/ plants11152058

Academic Editors: Milan S. Stankovic, Paula Baptista and Petronia Carillo

Received: 30 June 2022 Accepted: 1 August 2022 Published: 6 August 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Seasonal pHx changes are likely the consequence of steadied water and nutrient transport in the xylem that, in turn, is under environmental control and affected by winter stress. Under laboratory and field conditions, several abiotic and biotic stresses have been shown to cause an alkalinisation of the xylem or apoplastic sap, including drought [18–20], salt stress [21,22], flooding [23,24], chilling temperatures [25] and fungal infection [25,26]. Additionally, it has been shown that environmental factors that influence transpiration rates, such as increased solar irradiation, vapour pressure deficit (VPD) and increased temperatures can influence the pHx [16,27,28]. Particularly in response to drought, alkalinisation of the xylem sap induces accumulation of the phytohormone abscisic acid (ABA) in the leaf apoplast triggering stomatal closure [18,20,29]. During winter, soil frosts can induce drought stress, which becomes even more severe at high elevations, where winter desiccation is a widespread phenomenon observed in woody plants [30]. Little is known about the effects of subalpine winter environmental conditions on the alkalinisation of the pHx.

The aims of the present study were (1) to determine seasonal changes of pHx in a deciduous angiosperm (apple tree, *Malus domestica*) and in two evergreen gymnosperm tree species (*Picea abies*, *Pinus cembra*) and to correlate pHx values with seasonal temperature changes, (2) to compare seasonal changes of pHx in the same species under contrasting environmental conditions and growing sites (greenhouse versus field conditions in apple trees in order to advance tree phenology, low versus high elevation in evergreen gymnosperms) and (3) to assess the short-term dynamic of pHx in relation to developmental phenology during flushing of apple trees in spring.

#### **2. Results**

The pHx of apple trees was found to change significantly with season, with a mean value of pH 5.0 ± 0.03 in spring and pH 5.6 ± 0.05 in winter (Figure 1). The range between the minimum (pH 4.7) and maximum pH value (pH 5.9) recorded during the two measurement years was 1.2 pH units. In winter, the xylem sap was more alkaline. Acidification of the xylem sap occurred during regrowth in spring and was strongly linked to the developmental stage of buds during flushing. While in the early developmental stages, "green tip" and "half inch green", the pH remained unchanged; in later stages, beginning with a "tight cluster", the pH of the xylem sap decreased significantly (*p* < 0.05). In all later developmental stages, pHx values were more acidic and significantly different from the early developmental stages.

**Figure 1.** Seasonal change of pHx of twigs of *M. domestica* as related to the developmental stage of buds during flushing in spring. The grey boxes indicate the minimum and maximum for the respective season. Data (mean values ± SE, replicate numbers are given in the figure) are from December

to May from two measurement years and were obtained on potted trees and orchard trees at field sites. Different letters indicate significant differences (tested by ANOVA and the Tamhane test at *p* < 0.05). The pictures illustrate the phenology or stage of bud development.

The onset of or release from soil frost had no effect on pHx in apple trees. The seasonally occurring temperature fluctuations (mean values of soil and air temperatures 7 days before sampling) did not have any immediate effect on pHx of apple trees (see Table 1).

**Table 1.** Effect of air and soil temperatures (mean values from 7 days before sampling) on pHx (mean values, *n* = 5). The Pearson correlation coefficient (r) is given for calculations obtained with data from eight different sampling dates.

<sup>1</sup> measured on twig samples in an apple orchard from October 2010 till May 2011; <sup>2</sup> measured on twigs sampled at 600 m a.s.l. from October 2008 till April 2009; <sup>3</sup> measured on twigs sampled at 1950 m a.s.l. from October 2008 till April 2009.

Similar seasonal pHx changes were obtained for potted apple trees exposed to two contrasting temperature conditions (Figure 2). The range between mean values measured in spring (pH 4.4 ± 0.02) and mean values measured in winter (pH 5.7 ± 0.1) was 1.3 pH units. In potted trees, pHx values in spring were slightly more acidic than in trees from the field site. In the trees exposed to greenhouse conditions xylem sap acidification and flushing occurred earlier. In the greenhouse, apple trees were already at the developmental stage of "first bloom" in the middle of March. In the field, potted trees started to bloom nearly one month later.

**Figure 2.** Seasonal change of pHx of twigs of potted *M. domestica* trees grown in an orchard (closed circles) or under greenhouse conditions (10/25 ◦C; open circles). Data are mean values ± SE (*n* = 4). Trees were transferred into the greenhouse on 15 January 2010.

In the investigated evergreen gymnosperms, a significant seasonal change of pHx was also found. In *P. abies* at 600 m a.s.l. alkalinisation in winter was more pronounced than in apple with mean values of pH 5.4 ± 0.08 in spring and pH 6.9 ± 0.13 in winter (Figure 3). The range between the minimum (pH 5.2) and maximum pH value (pH 7.3) measured was 2.1 pH units. In *P. cembra* at 600 m a.s.l., mean values of pH 6.1 ± 0.3 in spring and pH 6.8 ± 0.07 in winter were found (data not shown).

**Figure 3.** Seasonal change of pHx of twigs of *P. abies* sampled in the Botanical Garden in Innsbruck at 600 m a.s.l. Different letters indicate significant differences between mean values (±SE; *n* = 5) tested by ANOVA and the Bonferroni test at *p* < 0.05.

Soil and air temperatures at the two contrasting field sites clearly differed (Figure 4). Particularly at the high elevation site (1950 m a.s.l.), there was a prolonged soil frost period during winter. Winter xylem sap alkalinisation of *P. abies* and *P. cembra* from lowland (600 m a.s.l.) and from the timberline ecotone sites (1950 m a.s.l.) are compared in Figure 5. In winter, xylem sap was more alkaline at the growing site at 1950 m independent of species (*p* < 0.001). In both evergreen gymnosperms, a significant alkalinisation was observed during winter at both elevations, with a stronger alkalinisation at the high elevation growing site (*p* < 0.001). pHx values in the gymnosperms showed a negative correlation with the seasonally occurring air and soil temperatures (*P. abies* at 600 m and soil temperature r = −0.901, *p* < 0.05 see Table 1). In contrast to apple trees, the xylem sap of the evergreen gymnosperms acidified before bud break in spring.

**Figure 4.** Temperature conditions in the Botanical Garden at 600 m a.s.l. (**a**) and at the field site at the timberline at 1950 m a.s.l. on Mt. Patscherkofel (**b**). The solid line shows the daily mean soil temperature at 5 cm depth, the grey area marks the daily minimum and maximum air temperatures measured at a height of 2 m.

**Figure 5.** Seasonal change of pHx of (**a**) *P. abies* and (**b**) *P. cembra* twigs sampled at 600 m a.s.l. in comparison to samples from the alpine timberline at 1950 m. White boxes show values obtained during the growing period (*n* = 9), grey boxes show winter values (*n* = 25). The boxplots present the median and the 10th, 25th, 75th and 90th percentiles; outliers are shown as dots. For both species at both growing sites winter values were significantly different (Bonferroni at *p* < 0.05) from values obtained during the growing period. Winter alkalinisation was significantly higher (*t*-test at *p* < 0.001) at 1950 m for both species.

#### **3. Discussion**

In the three investigated species, the pHx showed a pronounced seasonal variation. While in spring and autumn, pHx was acidic, a significant alkalinisation took place during winter. The measured pHx values were in the range of other published values for pHx [1,19]. The alkalinisation measured over winter corroborates earlier observations for a number of temperate deciduous angiosperm trees (*Acer pseudoplatanus* pH 5.4–6.9 [3]; *Actinidia chinensis* pH 5.3–6.2 [4]; *Betula pendula* pH 5.7–7.5 [5]; *Populus* × *canadensis "robusta"* pH 5.4–7.5 [6]; *Fagus sylvatica* pH 4.8–6.7 [31]; *Robinia pseudoacacia* pH 5.2–6.0 [8]).

In spring and autumn pHx of the gymnosperms *P. cembra* and *P. abies* was acidic. Acidification was also reported for other gymnosperms (pH 5.3 in *Abies koreana* [32]; pH 5.6 in *Pinus taeda* [33]). The seasonal dynamic of pHx in gymnosperms has recently been studied [17]; as for deciduous angiosperms, we found a significant winter alkalinisation in evergreen gymnosperms. The seasonal pHx amplitude in the two studied evergreen gymnosperms was higher than in apple trees, but still in the range of the reported maximum seasonal amplitude of 2.3 pH units as reported for *Fagus sylvatica* [7,31]. Our results allow us to compare between the studied species but not between the functional groups of trees. Further studies with a higher number of conifer and deciduous tree species are needed to compare between functional groups of trees at different environments.

In apple trees, xylem sap acidification in spring was linked to the bud developmental stage "tight cluster" and there was no significant correlation between xylem sap acidification and the seasonally occurring soil or air temperatures. A close relationship between developmental stages during flushing and pHx values was found in trees grown under field conditions and individuals exposed to greenhouse conditions. At the onset of bud break, when the release from winter dormancy first becomes visible, no significant changes of pHx values were found. In the developmental stage of "tight cluster" when foliation had proceeded to such an extent that a significant transpiration and consequent ascent of water and nutrients in the xylem were possible, pHx began to decrease. Concomitant measurements of transpiration and xylem conductivity reveal that in the developmental stage of "tight cluster" the developing leaves are already transpiring; however, xylem hydraulic conductance is not yet fully restored [34,35].

For apple trees, significant seasonal variations in the quantitative mineral and amino acid composition of the xylem sap are reported [36,37]. Moreover, for *P. abies*, the mineral composition of the xylem sap and its seasonal variation are reported [38,39]. In the investigated evergreen gymnosperms xylem sap acidification in spring occurred before bud break and therefore before the new needles of the current season emerged. We assume that in evergreen species the absence of water movement in the xylem in winter is directly related to winter stress and low environmental temperatures; therefore, a correlation between pHx values and temperatures can be found. In deciduous tree species, transpiration and water movement in the xylem can start only after foliation in spring and therefore pHx values might not be directly related to the environmental temperatures in spring. Tree phenology is linked to occurring environmental temperatures and therefore climate change can lead to phenological shifts [40]. These shifts in tree phenology will directly influence the timing of the seasonal variations in pHx for apple trees. In the two studied conifers at the alpine timberline, a correlation between pHx and the seasonal occurring soil temperatures was found. Increasing temperatures due to climate change will therefore also affect the seasonal changes of pHx in timberline conifers.

The underlying mechanisms for the seasonal variation of pHx are not fully understood. All factors affecting proton pumping activity might be involved [11]. Fromard et al. [8] showed that in *Robinia pseudoacacia* the plasma membrane H+-ATPase of the vessel-associated cells is responsible for the control of pHx, and its activity changes with season. In winter the activity of the H+-ATPase was low and therefore xylem sap was more alkaline. Seasonal changes in the amount of the H+-ATPase in cambial and expanding xylem cells are reported for twigs of two *Populus* species [12]. Furthermore, the xylem sap composition has been shown to influence the pHx due to the buffering capacity of distinct components [14,15,41,42]. In the xylem sap of beech roots, a significant correlation between low pH values and high concentrations of Ca2+, Mg2+ and malate was found [13]. The author concluded that the increase in malic acid, which forms complexes with cations, is responsible for the springtime acidification of the xylem sap in beech roots.

In the case of the winter acidification of the xylem sap of *Juglans regia* [10], the pHx values were not directly related to the activity of the H+-ATPase of vessel-associated cells, but depended on the seasonal variation of the sugar content in the xylem sap. Protoncoupled active sugar transport mechanisms were shown to be responsible for this [10,43,44]. There is a strong relationship between pHx and the concentration of dissolved CO2 in the xylem sap. Therefore, pHx regulates the amount of CO2 transported through the xylem and, furthermore, the exchange of CO2 between the different stem tissues [1,45,46]. Knowledge about the species-specific xylem sap pH might be useful for the calculation of CO2 budgets of individual trees.

The comparison of pHx values of the same species from contrasting elevations clearly revealed that at 1950 m a.s.l. winter alkalinisation was more pronounced than at 600 m in both investigated gymnosperms. The amount of mean winter alkalinisation in *P. abies* was 1.5 pH units at 600 m compared to 2.0 pH units at 1950 m. Similarly, in *P. cembra* we measured a mean winter alkalinisation of 0.6 pH units at 600 m compared to 1.5 pH units at 1950 m. In spring, pHx values were unaffected by elevation in the two gymnosperms. In both species at both growing sites, a significant correlation between pHx values and soil temperatures was found. At the timberline during winter the soil was frozen for a prolonged period causing winter drought stress that increased in severity with duration [30]. Timberline conifers are growing at the upper distribution boundaries of the respective tree species. Information on pHx and xylem sap composition may help to understand the survival mechanisms of the respective tree species at the timberline [17]. Under laboratory and field conditions, drought has been repeatedly identified as a factor in xylem sap alkalinisation, due to a reduced proton pumping activity, leading to the accumulation of abscisic acid (ABA) in the leaf apoplast and inducing stomatal closure [11,18,27,28,47]. The alkalinisation of the xylem sap in response to drought stress is not a universal mechanism in higher plant species [32]. However, the response to the apoplastic alkalinisation elevated ABA concentration in the apoplast and induction of stomatal closure—is thought to be a universal mechanism in plants [32]. The winter alkalinisation of the xylem sap in evergreen gymnosperms might be part of a physiological stress response in evergreen trees, which keeps the stomata securely closed during winter. Nevertheless, irrespective of its physiological function, pHx can be nicely used to assess the current activation state of the xylem tissue in trees.

#### **4. Materials and Methods**

#### *4.1. Study Site and Plant Material*

Seasonal changes of the pHx were studied in a deciduous (apple tree, *Malus domestica* Borkh. cv. "Golden Delicious" growing on "M9" rootstock) and two evergreen tree species (*Picea abies* L. Karst. and *Pinus cembra* L.). Samples from apple trees were either taken from trees growing in an apple orchard in Tarsch, Italy (46◦36 N, 10◦53 E, 860 m a.s.l., 20 trees were used for sampling) or from 20 potted trees cultivated in the Botanical Garden, University of Innsbruck, Austria (47◦16 N, 11◦23 E, 600 m a.s.l.). To break dormancy ahead of field sites and to advance tree phenology, four potted apple trees were transferred into a greenhouse (10/25 ◦C night/day) on 15 January 2010. pHx values of these trees growing in the greenhouse were then compared to potted trees kept under field conditions. Apple trees were from 5 to 6 years old, and 2-year-old and 3-year-old shoots were used for xylem sap extraction. Samples from the evergreen gymnosperms were taken from two different field sites with contrasting elevations: at 600 m a.s.l. from trees growing in the Botanical Garden of the University of Innsbruck, and from trees growing close to the timberline ecotone on Mt. Patscherkofel (1950 m a.s.l.; 47◦12 N, 11◦27 E). Xylem sap samples from *P. abies* and *P. cembra* (twigs with a diameter of 0.7–0.9 cm, detached from adult trees) were

taken from the end of October 2008 to the end of May 2009. Apple trees were sampled over two years from December 2009 to May 2011. Twig samples were always taken between 10.00 and 12.00 a.m.

#### *4.2. Temperature Conditions*

At the field sites, soil and air temperatures were recorded with sets of type T thermocouples. Soil temperatures were measured at a depth of 5, 10 and 20 cm, air temperatures were measured at a height of 2 m. Thermocouples were connected to a multiplexer (AM16/32B, Campbell Scientific, Logan, UT, USA) and temperatures were recorded at 7 min intervals with a CR10 data logger (Campbell Scientific, Logan, UT, USA). Temperatures of the potted apple trees in the greenhouse experiment and in the field were measured every 30 min with TidBit temperature data loggers (Onset, Pocasset, MA, USA).

#### *4.3. Xylem Sap Collection*

For xylem sap extraction, twigs with a diameter of 0.7–0.9 cm were used. Twig samples were detached from the trees and transported in plastic bags from the field sites to the laboratory. Samples were either prepared for immediate measurement or stored overnight in a cold room at +4 ◦C until measurements started. When twigs from the alpine timberline site were frozen at the time of sampling, they were stored overnight at −8 ◦C in commercial freezers until the beginning of measurements. Xylem sap was collected using a Scholander pressure bomb (Model 3115, Soil moisture Equipment Corp., Santa Barbara, CA, USA). Twigs were cut into pieces of a mean length of about 15 cm. The bark, including the phloem and the cambium layer, was peeled off to ensure that the expressed liquid came only from the xylem. The pieces of wood were then sealed into the lid of the pressure chamber. After closure of the pressure chamber, the pressure was gradually increased to a maximum pressure of between 0.2 and 2.5 MPa. Xylem sap leaking out at the cut surface of the wood was then collected with glass capillaries until the volume was sufficient for pH measurements. pHx was measured immediately after xylem sap collection. This sampling methodology may affect xylem sap composition [48]. To minimise this aspect, we applied only the minimum pressure necessary for xylem sap collection and worked with pH microsensors that allowed measurements with small quantities of xylem sap (see following paragraph).

#### *4.4. pH Measurements*

In apple twigs the pHx of the expressed xylem sap was measured with two different pH sensor types.  The pH-1 micro was connected to a PC for digital registration of the pH values. Because of the small dimension of the sensor (sensor tip < 150 μm), measurements were possible even within sap volumes < 0.02 mL. The second sensor type employed to measure apple xylem sap pH was a PHR-146 micro electrode (Lazar Research Laboratories, Los Angeles, CA, USA). For measurements with this sensor type the extracted xylem sap was also collected within glass capillaries, but then transferred to microtiter plates and measured in the microtiter plate. One droplet of xylem sap (0.05 mL) was sufficient for pH measurement. The pH sensor was connected to a digital pH meter (Jenco Model 60, Jenco Instruments, San Diego, CA, USA) and pH values could be read from the display. Calibration of the pH sensors was performed with certified buffer solutions (Merck, Darmstadt, Germany) according to the user instructions of the respective pH sensor. Additionally, all recordings were cross checked with non bleeding pH-indicator strips (Merck, Darmstadt, Germany). Measurements of the xylem sap of gymnosperms were taken using a micro combination pH electrode (Amani-1000, Innovative Instruments, Tampa, FL, USA). This sensor has a tip diameter of 1 mm and was also injected into the glass capillaries (inner diameter 1.6 mm) containing the collected xylem sap. The sensor was connected to a data logger (CR10X, Campbell Scientific, Logan, UT, USA) for digital registration of the values. Between 4 and 5 twig samples were used for the pHx measurements.

#### *4.5. Tree Phenology*

Classification of floral bud phenology in apple trees was conducted according to the descriptions given in the BBCH scale [49]. In the text the following terminology is used: "dormancy" (BBCH 00), "green tip" (BBCH 53), "half inch green" (BBCH 54), "tight cluster" where leaves are unfolding (BBCH 56), "first bloom" (BBCH 60) and "post bloom" (BBCH 69). In the two studied conifer species the timing of bud break in spring was observed.

#### *4.6. Statistical Data Analysis*

pHx values are given as mean ± standard error of the mean (SE). After values passed the Kolmogorov–Smirnoff and the Levene tests, significant differences between mean values were tested using one-way ANOVA and the Bonferroni post hoc test (*p* < 0.05). If homogeneity of variances was not established, significant differences between mean values were tested with one-way ANOVA and the Tamhane test. pHx values of gymnosperms from the two contrasting elevations were compared using Student's *t*-test. Correlations between mean pHx values and the seasonally occurring air and soil temperatures were calculated using Pearson's correlation analysis. All analyses are carried out using PASW Statistics 18 (formerly SPSS, IBM Corporation, New York, NY, USA). Replicate numbers are given in the figures or in the figure legends.

**Author Contributions:** Conceptualization, G.N. and M.P.; methodology, M.P. and E.L.; formal analysis, M.P. and E.L.; investigation, M.P. and E.L.; resources, G.N.; writing—original draft preparation, M.P.; writing—review and editing, G.N. and M.P.; funding acquisition, M.P. and G.N. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was financially supported by 'Forschungsförderungsmittel der Südtiroler Landesregierung' [PJ 81817] and a grant from 'Verein zur Förderung der wissenschaftlichen Ausbildung und Tätigkeit von Südtirolern an der Landesuniversität Innsbruck'. The authors thank the Department of Innovation, Research, University and Museums of the Autonomous Province of Bozen/Bolzano for covering the Open Access publication costs.

**Acknowledgments:** We thank M. Abler from 'Südtiroler Beratungsring für Obst- und Weinbau' and O. Buchner for assistance at the field sites. Further, we thank the Patscherkofel Bergbahnen for free transportation to Mt. Patscherkofel.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

### *Article* **Genome-Wide Identification and Spatial Expression Analysis of Histone Modification Gene Families in the Rubber Dandelion** *Taraxacum kok-saghyz*

**Francesco Panara 1, Carlo Fasano 1, Loredana Lopez 1, Andrea Porceddu 2, Paolo Facella 1, Elio Fantini 1, Loretta Daddiego <sup>1</sup> and Giorgio Perrella 3,\***

**Abstract:** *Taraxacum kok-saghyz* (*Tks*), also known as the Russian dandelion, is a recognized alternative source of natural rubber quite comparable, for quality and use, to the one obtained from the socalled rubber tree, *Hevea brasiliensis*. In addition to that, *Tks* roots produce several other compounds, including inulin, whose use in pharmaceutical and dietary products is quite extensive. Histonemodifying genes (HMGs) catalyze a series of post-translational modifications that affect chromatin organization and conformation, which, in turn, regulate many downstream processes, including gene expression. In this study, we present the first analysis of HMGs in *Tks*. Altogether, we identified 154 putative *Tks* homologs: 60 HMTs, 34 HDMs, 42 HATs, and 18 HDACs. Interestingly, whilst most of the classes showed similar numbers in other plant species, including *M. truncatula* and *A. thaliana*, HATs and HMT-PRMTs were indeed more abundant in *Tks*. Composition and structure analysis of *Tks* HMG proteins showed, for some classes, the presence of novel domains, suggesting a divergence from the canonical HMG model. The analysis of publicly available transcriptome datasets, combined with spatial expression of different developmental tissues, allowed us to identify several HMGs with a putative role in metabolite biosynthesis. Overall, our work describes HMG genomic organization and sets the premises for the functional characterization of epigenetic modifications in rubber-producing plants.

**Keywords:** histone modification; gene expression; *Taraxacum kok-saghyz*; natural rubber

#### **1. Introduction**

The genus *Taraxacum Wigg*. (dandelion) is included in the large family of Asteraceae (Compositae), which counts more than 2800 species [1]. Most of these species of dandelions live in the temperate zones of the northern hemisphere, although they are native to Eurasia [2]. Among all the dandelions species, the Russian dandelion (*Taraxacum kok-saghyz*, *Tks*) can produce from its roots natural rubber (NR) of excellent quality, quite comparable to that obtained from *Hevea brasiliensis*, the so-called rubber tree [3]. In *Tks*, the length of the rubber polymer is greater than that of Hevea, with a molecular weight of approximately 2180 kDa, but its productivity per hectare is lower [4]. Nevertheless, *Tks* mature roots can contain up to 5% of rubber and 20% of inulin in dry weight in wild plants and, therefore, this dandelion is considered as a valid alternative natural rubber source. *Tks* can extensively grow both in cold and temperate areas, presents a short life cycle, and can be easily harvested [5]. Moreover, *Tks* is a species suitable for genetic manipulation and it could be used as a model plant to investigate the molecular pathways that regulate NR biosynthesis [6,7].

NR is a polymer composed of cis-1,4-polyisoprene and it represents a very important raw material, used to produce more than 50,000 both industrial and medical products [8].

**Citation:** Panara, F.; Fasano, C.; Lopez, L.; Porceddu, A.; Facella, P.; Fantini, E.; Daddiego, L.; Perrella, G. Genome-Wide Identification and Spatial Expression Analysis of Histone Modification Gene Families in the Rubber Dandelion *Taraxacum kok-saghyz*. *Plants* **2022**, *11*, 2077. https://doi.org/10.3390/ plants11162077

Academic Editor: Pavel Kerchev

Received: 1 July 2022 Accepted: 4 August 2022 Published: 9 August 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The *Tks* genome has been recently published, showing that this dandelion (a diploid species 2 *n* = 16) has a genome size of 1.29 Gb, comprising 46,731 predicted protein-coding genes [9].

*Tks* presents tissues specialized in the production of NR, called laticifers [3]. Laticifers contain latex in specific structures, the rubber particles, that are spherical organelles surrounded by a lipid monolayer membrane [8]. Several proteins are bound to the membrane of the rubber particles and some of them play an important role in rubber biosynthesis. In the roots of Taraxacum spp., there are also present in significant quantities the carbohydrate inulin and other metabolites that are relevant in various applications [2,10–13].

In plants, chromatin is comprised in nucleosomes that are constituted by two copies of each histone protein (H2A, H2B, H3, and H4) wrapped around ~146 base pairs of DNA [14]. Together with canonical histone proteins, nucleosomes can also contain histone variants that are usually integrated in specific chromatin regions and in response to external stimuli [15]. Histones are composed of N-terminal tails that are enriched in basic amino acids that are subjected to different post-translational modifications (PTMs) [16]. The simultaneous presence of different modifications defines the so-called histone code that is superimposed over the genetic code to regulate most cellular mechanisms [17]. Histone H3 and H4 tails can indeed be methylated, acetylated, and phosphorylated at different levels (mono-di-and tri). Interestingly, different lysine and arginine residues can be subjected to different modifications that ultimately work as positive or negative regulators of gene expression [18]. Histone methylation is catalyzed by histone methyltransferases (HMTs) and can occur primarily at H3 Lys4 (K4), Lys9 (K9), Lys27 (K27), and Lys36 (K36) [19]. HMTs comprise mostly of the SET domain group (SDG) enzymes that are divided in seven different classes: (1) E(Z) (enhancer of zeste) homologs, (2) ASH1 (absent, small, or homeotic discs 1) groups (ASH1 homologs (ASHH) and ASH1-related proteins (ASHR)), (3) TRX (Trithorax) groups (TRX homologs and TRX-related proteins), (4) SET and PhD domain, (5) SU(VAR)3–9 groups (together with SU(VAR)3–9 homologs (SUVH) and SU(VAR)3–9 related proteins (SUVR)), (6) truncated SET domain class and (7) Rubisco large (RBLSMT) and small (RBSSMT) subunit methyltransferases SET-related class [20–22]. In addition to that, plants also do present a second HMT family that includes the protein arginine methyl-transferases that contain the PRMA domain (PRMTs) [23].

While histone methylation is catalyzed by HMTs', histone demethylases (HDMs) are responsible for its removal [24]. In plants, there are two main classes of HDM enzymes that are distinguished through their mechanism of action. Indeed, the KDM1/LSD1-likes (HDMAs) operate through amine oxidation to remove methylation marks, while JmjC domain-containing (JMJs) demethylases use hydroxylation [25–27]. At the same time, cofactors and the substrates change based on the type of HDMs. HDMAs only demethylate residues that have been mono- or dimethylated by a flavin adenine dinucleotide (FAD) dependent reaction. Instead, JMJs demethylate lysines independently on the type of methylation and use Fe(II) and α-ketoglutarate (αKG) as cofactors [28].

Histone acetyltransferases (HATs) and deacetylases (HDACs) catalyze the addition and the removal of the acetyl-CoA groups over the histone tails, respectively [29]. This modification is usually correlated with changes in gene expression as well as chromosome condensation. In Arabidopsis, HATs are divided in four groups: GNAT and MYST families (HAGs and HAMs, respectively), CBP (HACs), and the TAFII250 (HAFs) [30]. HAGs are then classified in three main classes: GCN5 or GCN5-likes, ELP3, and HAT1. MYSTs, instead, contain mostly HAM members [31]. HDACs in plants are grouped in three main groups: RPD3/HDA1 superfamily, similar to the yeast large RPD3 complex (HDAs), the sirtuins (SIR2/SRTs), and the HD2/HDT enzymes, whose class appears to be specific in plants [32,33]. Both HATs and HDACs have been largely characterized in plants. Indeed, their function is involved in different developmental processes and transitions, including germination, cell differentiation, and leaf and floral organogenesis [34–37]. Furthermore, their role has been linked with changes in gene expression upon different abiotic and biotic stress conditions, such as salinity, light, temperature, and immuno-responses [38–42].

Here, we report for the first time an extensive in silico analysis and identification of the histone-modifying genes (HMGs) in *Tks*. Using publicly available genome and transcriptome data, we analyzed the gene and protein structure of the identified HMGs and monitored their expression in various *Tks* tissues, including latex. In addition, we measured transcript levels in *Tks* young and adult leaves and roots and identified those per each class of enzymes that appeared to be abundant in relevant tissues.

#### **2. Results**

#### *2.1. HMG Genes Identification*

In the *Tks* genome [9]*,* we identified a total of 154 HMGs: 60 HMTs, 34 HDMs, 42 HATs, and 18 HDACs. We compared the total number with a close relative of *Tks*, *Lactuca sativa* (Ls) and with the reference plants *Arabidopsis thaliana* (At) [43] and *Medicago truncatula* (Mt) [44]. While some classes were similar in numbers, others showed divergent expansion among different species. PRMTs were indeed more abundant in *Tks*, indicating an expansion of this class of proteins similarly to *Litchi chinensis* [45]. Altogether, HAGs were more abundant in both *Tks* and *Ls* compared to *At*, confirming the possibility for this group to expand divergently in different taxa (Table 1). The complete list of *Tks* HMGs is reported in Table S1.

**Table 1.** Number of genes identified in the different classes of HMGs in *Taraxacum kok-saghyz* and in three reference species: *Lactuca sativa*, *Arabidopsis thaliana*, and *Medicago truncatula*.

Gene duplication events contributed to the high number of HMGs in the *Tks* genome. Synteny analysis on the *Tks* genome assembly was not effective for the identification of segmentally duplicated HMGs, while five couples of tandemly duplicated genes were identified: *Tks*HDA5/*Tks*HDA6, *Tks*HAG30/*Tks*HAG31, *Tks*HAG36/*Tks*HAG37, *Tks*JMJ11/*Tks*JMJ12, and *Tks*SDG44/GWHPAAAA005121. In the last pair, GWHPAAAA005121 was not classified as an HMG due to the absence of a representative domain. We identified 32 additional couples of putative duplicated HMGs based on phylogenetic analysis. The Ka/Ks ratio was calculated for all the gene pairs in order to estimate the occurring evolutionary dynamics (Table S3). Most of the couples showed Ka/Ks <1, suggesting purifying or stabilizing selection. Only the *Tks*HAG13/14 pair with a Ka/Ks ratio of 1.54 and *Tks*PRMT7/8 with 1.05 showed values compatible with positive and neutral selection.

#### 2.1.1. HMTs

Among HMTs, 46 were identified belonging to SDG and 14 to the PRMT group. SDGs are similar in number to *At* and the same amount in *Tks* and *Ls* (Table 1). SDGs are divided in seven classes [22,46]. Four *Tks*SDGs cluster with Class I, E(Z)-like SDGs (Figure 1A) and show the expected SANT-CXC-SET domains except for *Tks*SDG8 lacking the SANT *domain* (Figure 2).

Five clusters with Class II, ASH1-like. *Tks*SDG29, *Tks*SDG16, and *Tks*SDG19 have the same domain composition of corresponding *At* proteins. *Tks*SDG15 encodes for a short (52aa) and a probably aberrant protein where only the SET domain was identified. *Tks*SDG37 lacks the PostSET domain compared to AtSDG26. Five proteins cluster with Class III, TRX-like. *Tks*SDG41 differs from AtSDG14 as the first PHD module is substituted by a SAND domain that is often associated to PHDs and could contribute to their function. *Tks*SDG45 showed an additional TUDOR domain at the N-terminal. The TUDOR domain was previously identified in *D. melanogaster* TUDOR proteins and its function is unknown. Several human JMJ/KDM4 proteins harbour a TUDOR domain at the C-terminal [44]. *Tks*SDG4 clusters with *At*SDG25 but differs in domain composition and sequence length: GYF-SET and 1057aa the former, and SET-PostSET 1424aa the latter. This was previously observed in *S. lycopersicum* (*Sl*), where *At*SDG25 and *Sl*SDG20 showed similar differences [45]. Two Class IV SDGs exist in both *Tks* and *At*. Six *Tks* proteins cluster with class V subclass I and six with class V subclass II with a similar domain composition to *At* proteins. Eighteen SDGs can be classified as belonging to class VI/VII and show an interrupted SET domain. *Tks*SDG40 is a short protein of 158aa, probably aberrant, and the presence of the SET domain was not confirmed by SMART analysis.

*Tks*SDG5 and *Tks*SDG25 encode for putative long proteins with an additional/unexpected domain composition, probably resulting from exon gain or, more likely, a defect in genome assembly/annotation.

As shown in Figure S1, three PRMTs form a separate cluster: *Tks*PRMT2, *Tks*PRMT7, and *Tks*PRMT8 and contains two PRMT domains (Figure 2). The PRMT5 domain was not confirmed for *Tks*PRMT9 and *Tks*PRMT11 by SMART analysis. *Tks*PRMT6 shows four C2H2 modules at the N-terminal similarly to HMGs belonging to other groups such as *At*SDG6/*Tks*SDG24 (Figure 2).

**Figure 1.** Phylogenetic tree of SDG (**A**) and JMJ (**B**) proteins of *Taraxacum kok-saghyz* and *Arabidopsis thaliana* (in colour according to their respective classes). The numbers near the tree branches represent bootstrap values.

**Figure 2.** Domain composition and intron-exon structure of *Taraxacum kok-saghyz* HMTs and HDMs.

#### 2.1.2. HDMs

In *Tks*, we identified 27 proteins containing the JmjC domain (PF02373) (Table 1). Despite that nine harbour the JmjC domain alone and could be classified as JMJ-only, only *Tks*JMJ25 and *Tks*JMJ27 clustered with *At* JMJ-only proteins (Figures 1B and 2). The other

seven are closer to other JMJ groups but underwent the loss of representative domains. In addition, only four out of fourteen proteins show the presence of the RING domain typical of the KDM3 group (Figure 2). KDM4 JMJs present N-terminal JmjN (PF02375) and JmjC, and C-terminal C5HC2 (PF02928) (subgroup I) or C2H2 (PF00096) (subgroup II). In *Tks*, two proteins belong to subgroup I and three to subgroup II. Four proteins cluster with KDM5 JMJs and show the expected domain composition (Figures 1B and 2).

In addition, among *Tks*JMJs, several proteins revealed additional C-terminal domains (*Tks*JMJ19, *Tks*JMJ3, *Tks*JMJ16, and *Tks*JMJ26) or exceptionally long introns (*Tks*JMJ7) that likely derive from errors in genome annotation (Figure 2). In the JMJD6 group, *Tks*JMJ15 has the expected F-Box domain while *Tks*JMJ5, although clustering with AtJMJ21, retains only the JmjC domain (Figures 1B and 2).

Seven HDMAs were identified by our analysis (Table 1). *Tks*HDMA1, *Tks*HDMA7, and *Tks*HDMA4 are grouped with *At*HDMA4, 2, and 1, respectively (Figure S2). The abovementioned proteins and *Tks*HDMA6 show the presence of SWIRM-NADB8-AOD domains (Figure 2). *Tks*HDMA2, *Tks*HDMA3, and *Tks*HDMA5 form a subgroup harbouring a SANT domain (Figures S2 and 2).

#### 2.1.3. HATs

We identified 42 HATs in *Tks*. Most of them, 38, belong to the HAG group (Table 1). An expansion of HAGs compared to *At* was already observed in other species analysed so far such as *Medicago truncatula* [44], *Malus domestica* [26], *Solanum lycopersicum* [47], *Citrus sinensis* [48]*,* and *Vitis vinifera* [49].

GCN5-type HAGs, *Tks*HAG36, and *Tks*HAG37 have a C-terminal Bromodomain similar to *At*HAG1. Three additional proteins clustering with *At*HAG1 are shorter and probably underwent loss of the Bromodomain: *Tks*HAG8 193 aa, *Tks*HAG33 250 aa, and *Tks*HAG34 192 aa (Figures 3A and 4 and Table S1).

**Figure 3.** Phylogenetic tree of HAG (**A**) and HDA (**B**) proteins of *Taraxacum kok-saghyz* and *Arabidopsis thaliana* (in colour according to their respective classes). The numbers near the tree branches represent bootstrap values.

**Figure 4.** Domain composition and intron–exon structure of *Taraxacum kok-saghyz* HAGs.

*Tks*HAG16 clusters with *At*HAG2 and harbours the Hat1 domain. *Tks*HAG10 clusters with *At*HAG3 and harbours the ELP3 domain. Two additional proteins clustering with *At*HAG3, *Tks*HAG27, and *Tks*HAG6 show the AT1 domain alone (Figures 3A and 4). Furthermore, other AT1 domain-containing proteins can be observed.

A first group of seven proteins forms a definite cluster that is characterized by the presence of Jas and PHD domains except for two shorter proteins, probably derived from the others, that lost specific domains: *Tks*HAG21 and *Tks*HAG25. In *Tks*HAG32, the AT1 domain was not identified as below the threshold for SMART analysis (Figures 3A and 4).

Ten proteins are characterized by the presence of an additional FR47 domain at the C-terminal. Four of them, *Tks*HAG5, *Tks*HAG11, *Tks*HAG17, and *Tks*HAG19 form a definite cluster. Instead, a second cluster includes *Tks*HAG14, *Tks*HAG15, *Tks*HAG23, *Tks*HAG29, and two additional proteins: *Tks*HAG20 where no domain was identified by SMART, and *Tks*HAG24, where the FR47 domain is missing and there is an N-terminal AAK domain. Interestingly, *Tks*HAG13, containing both AT1-FR47, do not cluster with the other FR47 domain-containing proteins (Figures 3A and 4).

Two proteins, *Tks*HAG9 and *Tks*HAG26 show two C-terminal BRCT domains and form a cluster with *Tks*HAG30 and *Tks*HAG31 that are shorter and without a BRCT domain. *Tks*HAG12 harbours AT1 and C2 domains. All other *Tks*HAG proteins show the AT1 domain alone (Figures 3A and 4).

One MYST acetyltransferase was identified: *Tks*HAM1 with Chromo-C2H2-MYST-Syja\_N domain composition (Figure S3).

Two HACs identified in *Tks* have a TAZ-PHD-KAT11-ZZ-ZZ-TAZ domain composition, a third HAC, *Tks*HAC2, lost the N-terminal TAZ domain and clusters with *At*HAC2 (Figures S3 and S4). No proteins showing the HAF domain were identified by our analysis.

#### 2.1.4. HDACs

We identified 18 histone deacetylases (Table 1). Nine are HDAs, three belong to class I, RPD3-like, five to class II, HDAC1-like, and one to class III, HDAC11-like. In class I, *Tks*HDA9, clustering with *At*HDA19, is a long protein (1634 aa) showing six additional N-terminal TPR domains probably resulting from exon gain or a defect in genome assembly/annotation. In the class II, *Tks*HDA7 similarly to *Mt*HDA8 [44] harbours an N-terminal Znf\_RBZ domain (Figures 3B and S5).

Four proteins are SRTs containing the SIR2 (PF02146) domain. Interestingly, with the exception of *Tks*SRT4, all other SRTs have a double SIR2 domain (Figures S5 and S6).

Five *Tks*HDTs were identified by BLASTp analysis (Table 1). The most conserved are *Tks*HDT1 and *Tks*HDT2 if compared with *At*HDTs. *Tks*HDT1 showed a C2H2 domain. In *Tks*HDT5, C2H2 modules are six (Figures S5 and S7). In plants, HD2-Types are related to cis–trans isomerases found in insects and yeast [30]. Similarly, their function as HDACs remains debatable.

#### *2.2. In Silico Expression Analysis of HMG Genes in Different Organs and Developmental Stages*

To shed light on the spatial expression of the HMGs in *Tks*, we analyzed publicly available datasets [9], based on the following tissues/organs: leaves, stems, the main root, and the lateral root in both young and mature stages (eight samples); latex and three additional samples from reproductive organs: flower, peduncle, and seeds. Expression patterns of HMGs were evaluated for each gene as log2 fold change (Fc) compared to the average among samples. Hierarchical clustering was used to assess groups with similar expression patterns.

#### 2.2.1. HMTs

Among SDGs, a first group presents a peak of expression in young/mature leaves compared to other organs. Within this group, *TksSDG28* shows an additional peak in the latex. A second large group includes genes with higher expression in reproductive organs, mainly flowers and seeds. In group 3 we identified genes with higher expression in young organs and flowers. Group 4 shows higher expression in stem and roots with lower expression in the latex. The last group harbours the most expressed genes, showing an expression pattern primarily in mature roots (Figure 5).

*Tks*PRMTs exhibit less modulated expression patterns. Three groups can be identified: the first with higher expression in latex and reproductive organs; the second with a peak in flowers; and the third group with higher expression in stem and roots. In the first group, *TksPRMT7* is the one showing the most abundant peak of expression in latex (Figure 5).

#### 2.2.2. HDMs

Among *Tks*JMJs, a first group includes genes with higher expressions in reproductive organs. In addition to that, *TksJMJ22* presents a relevant peak also in the latex, indicating a strong modulation compared to whole roots. A second group includes genes with a higher expression in stem-root and flower, and little expression in latex and leaf. The third group shows a higher expression in stem-root but, differently from most JMJs, not in flower. This group includes JMJs with high levels of expression. The fourth group shows a relevant expression in stem and flower (Figure 5).

All HDMAs exhibit higher expression in stem, flower, and seeds, with two groups that can be distinguished by high and low expression in mature root and latex (Figure 5).

#### 2.2.3. HATs

Four groups were identified among *Tks*HAGs: a first one with peaks in the roots (mainly in mature roots); a second group predominantly in reproductive organs; a third with the same pattern with addition of leaves; and interestingly, the fourth group is made of highly expressed genes showing a peak in the latex similarly to the single HAM gene detected by our analysis. HACs are less expressed in leaves and latex (Figure 6).

**Figure 5.** Heat map of *Tks*HMTs and *Tks*HDMs in different organs and developmental stages. The main clusters are indicated by square brackets. Genes with higher expression are indicated in bold. Asterisks indicate the genes analyzed by qPCR.

**Figure 6.** Heat map of *Tks*HATs and *Tks*HDACs in different organs and developmental stages. The main clusters are indicated by square brackets. Genes with a higher expression are indicated in bold. Asterisks indicate the genes analysed by qPCR.

#### 2.2.4. HDACs

*Tks*HDAs can be divided into three groups: a first one with peaks of expression in leaves and peduncle; a second group with lower levels in the roots; and a third group more relevant in young and reproductive organs. *TksHDA8* shows a peculiar pattern with peaks in the main root and latex. *Tks*SRTs were present in reproductive organs with *TksSRT4* and *TksSRT2* having an additional peak in mature leaf and in the latex, respectively. HDTs, instead, show peaks of expression in young leaf, stem, lateral root, and in the mature stem with the exception of *TksHDT5*, which shows peaks in reproductive organs.

#### *2.3. HMGs Expression*

Using the absolute quantification method, we measured the expression of 13 different *Tks*HMGs belonging to four families in four different tissues: developing leaves (DL), fully developed leaves (FL), root tips (RT), and mature roots (R). Figure 7 shows the spatial gene expression pattern of four HMTs (*TksSDG3, TksSDG4, TksSDG24,* and *TksPRMT9*), five HDMs (*TksJMJ14, TksJMJ23, TksJMJ25, TksHDMA4*, and *TksHDMA6*), one HAT (*TksHAG10*), and three HDACs (*TksHDA7*, *TksHDA8*, and *TksHDA9*). Five genes (*TksHDMA4, TksHDA9, TksHDMA6, TksJMJ23*, and *TksHAG10*) have at least one tissue expression level higher than 1 femtogram (fg)/μL (Figure 7a). The most expressed gene is *TksHDMA4*, with an average tissue quantification of 6 fg/μL but no significative modulation between the different tissues. Among the remaining eight genes, *TksSDG3* and *TksJMJ25* expressions never exceed 0.1 fg/μL (Figure 7b). With an average tissue quantification of 0.05 and 0.01 fg/μL, respectively, they are the lowest expressed HMG genes. Overall, there is no difference of global HMG gene expression between the four different tissues analyzed, with an average expression of 1.38 and 1.29 fg/μL in fully developed leaves (FL) and the root (R), respectively. *TksJMJ25* and *TksHDA8* expression increased significantly, 319 and 174%, respectively, in FL in comparison to DL. *TksHAG10*, *TksSDG24*, and *TksPMRT9* are up-regulated in a less strong way. The remaining eight genes confirm an FL up-regulation with respect to DL; however, below a 50% change. *T*ks*SDG4*, *TksHDA8*, and *TksSDG3* expression increased clearly with a 276, 249, and 104%, respectively, in R compared to RT. *TksHDA7* and *TksHAG10* are up-regulated by 85 and 56%. As for leaves, all eight remaining genes are up-regulated in R tissue in comparison to RT with a percentage below 50%.

**Figure 7.** qPCR expression analysis of 13 HM genes in *Taraxacum kok-saghyz* developing leaves (DL), fully developed leaves (FL), root tips (RT), and mature roots (R). (**a**) Absolute quantification (mean ± SD fg/μL) of HM genes with an average expression above 1 fg/μL; (**b**) Absolute quantification (mean ± SD ag/μL) of HMG genes with an average expression below 1 fg/μL.

#### **3. Discussion**

#### *3.1. HMG Members in Taraxacum Kok-Saghyz*

Natural rubber represents a biopolymer of major importance due to its wide properties that cannot be found in synthetic material, including resilience, elasticity, and abrasion [50]. The demands to develop new sources of natural rubbers have quickly increased over the years due to the constant reductions of petroleum-based materials and the efforts towards renewables [3]. Currently, the main resource for commercial natural rubber, *Hevea brasiliensis*, has generated major concerns because of its high sensitivity towards pathogen infections [51]. Furthermore, latex can induce several allergic reactions due to their protein content [52]. As an alternative to *H. Brasiliensis*, two plant species are so far known to produce rubber in similar amounts and a high molecular weight: *Parthenium argentatum* Gray, also known as guayule, and *Taraxacum kok-saghyz*, the Russian dandelion [53]. *Tks* latex, for instance, is broadly investigated for non-medical applications such as tire production [54].

In *Tks*, rubber biosynthesis occurs through different pathways that are responsible for NR chain elongation, small rubber particles (SRPPs), and rubber elongation factors, respectively [7]. Therefore, to improve rubber quality and quantity, it becomes paramount to functionally characterize the genes involved in those pathways.

Histone modifications are well established to control gene expression [15]. Thus, they represent a promising strategy to manipulate different signalling pathways, unlike the gene-by-gene approach [55].

In this study, we presented for the first time a comprehensive analysis of the HMGs in *Tks*, identified by the presence of representative domains recognized through their hidden Markov profiles. HMGs sequences and/or composition were compared with the model species *A. thaliana* and *M. truncatula*, as well as the close relative *L. sativa*.

Overall, the identification recorded an increase in the number of HMG genes in *Tks* (154) compared to *A. thaliana* (102), which is similar to what is observed in *L. sativa* and *M. truncatula*. Due to limitations in the assembly, it was not possible to assess gene duplication.

A closer investigation revealed that such an increase did not occur on all the HMG classes, but it was indeed specific. Thus, while HDMs were similar in numbers (34 in *Tks* vs. 24 in Arabidopsis), HATs were predominantly higher in the former (42) than in the latter (12). This is not necessarily surprising, as an expansion of HAGs was already revealed in other crops when compared to Arabidopsis (Table 1) [44].

Among the HMTs family, analyses on *Tks* indicated the presence of 46 genes encoding for SDGs, a similar number to those of lettuce and Arabidopsis. Instead, the PRMTs were doubled up in *Tks* compared to the other species. Overall, the *Tks* SDGs showed a similar structure to those already annotated with SANT-CXC-SET domains except for *Tks*SGD8 that showed a shorter structure (Figure 2). SANT domains are essential to ensure histone substrates for the SGD enzymatic activity [56], while the CXC allows binding to RNA molecules as the SET represents the signature domain to cluster the SDGs [57,58]. *Tks*SDG5 appears to be much longer than the other SDGs, perhaps due to fusion between two different proteins (Figure 2). The presence of a bag6-A domain could implicate a divergent function among the other SDGs, including a possible interaction with heat shock proteins [59]. A functional study is therefore required to address this question.

Deeper investigation on the PRMTs indicated that *Tks*PRMT2, 7, and 8 contained two PRMT domains whilst the majority had only one. Structural analyses combined with domain deletion could retrieve more information on that regard.

*Tks* HDMs were classified in JMJs and HDMAs. Interestingly, many JMJ-only clustered with the KDM3 groups although they lacked the additional RING domain (Figure 5). However, they all presented a very extensive introns and exons structure that is typical of the KDM3s [60].

Over 90% of *Tks* HATs belonged to the HAG group. Similar to previous reports [48], the majority of the members contained only the AT1 protein domain (Figure 4). Whilst *Tks*HAG36 and 37 were identified as GCN5 likes due to the presence of the bromodomain, others included Jas and PhD domains or the FR47, or the BRCT domain in combination with

AT1. This appears to be quite divergent from the canonical HAGs structure, suggesting a novel class of HMG HATs in *Tks*. Indeed, Jas domains are present in JASMONATE ZIMdomain (JAZ) proteins that are repressors of jasmonate (JA) signalling [61]. FR47 domain resembles the C-terminal region of the Drosophila melanogaster hypothetical protein FR47. Interestingly, a member of the HATs in *Plasmodium falciparum* also contains that domain, suggesting a class of acetyltransferases targeting nonhistone proteins [62].

Our approach identified 18 HDACs in *Tks*. The three classes showed numbers comparable to Arabidopsis (Figure S5). Interestingly, *Tks*HDA9, which shows similarities with *At*HDA19, appears to have additional domains. *At*HDA19 is a major regulator of development and growth [41,63]. Whether *Tks*HDA9 has developed similar functions remains to be assessed.

#### *3.2. In Silico and Spatial Analysis of HMGs in Tks Tissues Identified Putative Regulators of NR Biosynthesis*

Histone modifications are well known to regulate important functions in plants [15]*,* but very little is known about their expression patterns and behaviour in *Tks*. To unveil new roles for the HMGs, we took advantage of publicly available datasets in *Tks* and measured their mRNA levels as fold-change normalized to the average amounts. Among the different tissues analyzed, latex represented the most interesting as it is directly connected to natural rubber production. Indeed, the biosynthesis of NR occurs in the latex of laticifers, organized in rubber particles [64]. Within the SDG class VI, *TksSDG28* presented a distinct expression in latex (Figure 5). Instead, *TksSDG3* and *4* showed an abundant level of transcripts in young leaves and roots. Interestingly, *TksSDG4* was more abundant in roots, whilst *TksSDG3* showed similar levels in both tissues. The Arabidopsis *Tks*SDG3 homolog, SWINGER, is part of a large protein complex that includes VRN2 (VERNALIZATION 2) and VIN3 (VERNALIZATION INSENSITIVE 3), together with FERTILIZATION INDEPENDENT ENDOSPERM (FIE) and CURLY LEAF (CLF). They are responsible for establishing FLC (FLOWERING LOCUS C) repression during vernalization [65]. Notably, the *Tks*SDG4 homolog in Arabidopsis is ATXR7 that instead promotes H3K4me3 on the FLC locus, hence its transcriptional activation [66]. A similar function is associated with *At*SDG8, a homolog of *Tks*SDG24 [67]. Whether their function is indeed conserved in *Tks* and is based on tissue specificity remains to be assessed.

Heatmaps for *Tks* PRMTs did not show major modulations among tissues, with the exception for *TksPRMT9* and *Tks*P*RMT11* that portrayed an elevated expression in young leaves (Figure 5). Real time qPCR analyses of some of the indicated genes revealed a tendency towards root samples. Interestingly, *TksSDG4* was expressed primarily in mature roots and *TksPMRT9* showed major peaks in fully developed leaves as well as roots, making them interesting targets for a functional approach (Figure 7b).

In silico expression of the histone demethylases revealed a predominance on the reproductive tissues (Figure 5). However, *TksJMJ14* showed a preferential expression in roots. Our qRT-PCR data showed instead that levels were different between root tips and adult roots, suggesting a modulation associated with plant development. Furthermore, it is important to highlight that, while mature roots contain laticifers already producing rubber, in the root tips the laticifers are still developing. Therefore, our qPCR analysis presents dual significance: first, it identifies putative HMGs involved in regulating the expression of genes related to rubber or laticifers production; second, the comparison between the expression in leaves and roots is pivotal to distinguish between overall growth and rubber-specific mechanisms. The analysis of other two JMJs *TksJMJ23*, a homolog of Arabidopsis ELF6 [67] and *TksJMJ25*, confirmed in silico results as none of them were modulated but instead they showed similar expression levels among the analyzed tissues (Figure 7a,b). A similar pattern was observed for the HDMAs, where only *TksHDMA4* turned out to be the most abundant among the analyzed genes. Interestingly, the heatmaps showed a major expression in the latex tissue.

Heatmaps for histone acetyltransferases indicated only a few HATs were expressed in root tissues, whilst the group 3 was majorly present in young and mature leaves (Figure 6). We focused our attention on *TksHAG10* that indeed showed to be mostly abundant in roots (Figure 7a). Interestingly, our qPCR data also showed a peak in fully developed leaves that was not observed in the datasets analysis. However, the plant material used in the datasets experiment might not fully correspond to the one obtained for expression analysis.

Among the histone deacetylases, *TksHDA9* was the most expressed within the family (Figure 7a). Interestingly, *TksHDA9* revealed a drop in expression between leaves and roots, while *TksHDA8* was instead more abundant in mature roots. It will be interesting to assess whether *TksHDA9* action is counteracted by histone acetyltransferases GCN5-likes.

#### **4. Materials and Methods**

#### *4.1. Plant Material*

A *Tks* plant belonging to the W6-35166 population obtained from USDA-ARS (Regional Plant Introduction Station, 295 CLARK HALL, WSU Pullman, WA 99164, Washington State University, USA) was grown and clonally propagated from root cuttings. Fully developed clones have been cultivated in trays filled with topsoil under controlled conditions (12/12 light cycle, 20 ◦C) for three weeks, then transplanted in 85 × 39 × 34 cm rectangular pots, 8 plants each. Potted plants were cultivated from January to November 2021 under near natural conditions for photoperiod, temperature, humidity, and exposure to atmospheric events, in ENEA Trisaia Research Center (40◦09 47.4" N 16◦38 00.1" E, Italy). Plants were supplemented with fertilisers, treated against pests and fungal pathogens, and shaded during summer. Four different tissues have been harvested and snap-frozen in liquid nitrogen: developing leaves, fully developed leaves, root tips, and mature roots.

#### *4.2. In Silico Identification and Analysis of HMG Loci*

To identify HMGs loci, the hidden Markov profiles of each gene family were used as a query input for the HMM software against the protein subject datasets of Taraxacum kok-saghyz and Lactuca sativa [9]. Since for the HDT family no PFAM domain is available, we used Arabidopsis HDTs as query input for BLASTp searches against the protein datasets. The identified loci were analysed with the SMART software to characterise the domain composition [68,69]. The cases pointing to two different domains for the same protein region were resolved based on the domain that was associated to the smaller E-value. The intron/exon structures were graphically represented with GSDS 2.0 [70].

#### *4.3. Phylogenetic Analysis*

The HM protein families identified in this work from Taraxacum kok-saghyz genome annotation database [9] along with HM proteins from Arabidopsis were aligned using ClustalW program in CLC Genomics Workbench version 9.5.2 (Qiagen, Hilden, Germany). A phylogenetic tree was then constructed using a neighbor-joining algorithm and Kimura protein substitution model. Reliability of the internal branching was obtained by a bootstrap test of 1000 replicates.

Collinearity analysis was performed using MCSCANX [71]. In brief, Taraxacum proteins were searched for sequence homology using BLASTp [72] according to McScanX default settings. Collinear blocks were identified by McSCANX with default settings and sequence duplication were classified with the duplicate gene classifier tool integrated in the McSCANX suite. Protein of duplicated HMG were aligned with ClustalO [73] and codon-based alignments were obtained with PRANK software [74]. The synonymous and non-synonymous substitution patterns were calculated with PAML [75] using YN00 procedure [76].

#### *4.4. Transcriptome Analysis*

RNA reads included in project PRJCA000437 were downloaded from the National Genomics Data Center (NGDC), part of the China National Center for Bioinformation (CNCB) and quality checked using Trimmomatic [77] according to Corchete et al. (2020) [78]. The reads were then aligned to the *Tks* indexed genome using STAR [79] with default settings. HTSeq-count script [80] was used to determine raw counts for each feature and experiment. The raw count matrix was normalised using edgeR and TMM method [78,81].

#### *4.5. RNA Extraction, cDNA Synthesis, and Gene Expression Analysis*

Total RNA from the above-mentioned plant material was extracted using innuPREP Plant RNA Kit (Analytik, Jena, Germany) according to the manufacturer's protocol. RNA quality and concentration were estimated by Nanodrop™ 1000 Spectrophotometer (Thermo Fisher Scientific Waltham, MA, USA). A total of 1 μg of total RNA was used for cDNA synthesis using SuperScript IV kit (Thermofisher Waltham, MA, USA) with oligo(dT)20 as primers in a 20 μL final volume. Primer pairs for SYBR-Green based qPCR have been designed for each candidate gene using Geneious Prime 2022.1.1 (https://www.geneious.com). Each primer has been subjected to off-target analysis using local BLAST-2.11.0+ [72,82,83] against the Taraxacum kok-saghyz genome annotation [9] downloaded from https://ngdc.cncb.ac.cn/gwh/. Pairs with both forward and reverse primers matching an off-target transcript were rejected. Primer sequences are reported in Supplementary Table S2. qPCR was performed using an ABI Prism 7900HT instrument (Applied Biosystems, Waltham, MA, USA) and Platinum SYBR Green qPCR SuperMix-UDG with ROX (Thermofisher, Waltham, MA, USA) following manufacturer's instructions. Reactions were performed in two technical replicates on three biological replicates. The following cycling conditions were used for quantitative PCR: 5 min at 95 ◦C followed by 45 cycles at 95 ◦C for 15 s and at 58 ◦C for 60 s. Melting curve analysis from 60 to 90 ◦C was performed to monitor the specificity of the amplification. mRNA levels of each genes analysed were calculated using absolute quantification based on the standard curve method [84] as reported by D'Amelia et al., (2014) [85] expect for standard curves obtained by purified and quantified conventional PCR normalised to a concentration of 33 pg/μL and 10-fold serial dilutions (ranging from 102 to 108). Expression data were analysed using Tukey's pairwise test.

**Supplementary Materials:** The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/plants11162077/s1, Figure S1: Phylogenetic tree of PRMT proteins of Taraxacum kok-saghyz and Arabidopsis (in bold). The numbers near the tree branches. Figure S2: Phylogenetic tree of HDMA proteins of Taraxacum kok-saghyz and Arabidopsis (in bold). The numbers near the tree branches represent bootstrap values. Figure S3: Domain composition and intron–exon structure of Taraxacum kok-saghyz HAMs and HACs. Figure S4: Phylogenetic tree of HAC proteins of Taraxacum kok-saghyz and Arabidopsis (in bold). The numbers near the tree branches represent bootstrap values. Figure S5: Domain composition and intron–exon structure of Taraxacum kok-saghyz HDAs, SRTs, and HDTs. Figure S6: Phylogenetic tree of SRT proteins of Taraxacum kok-saghyz and Arabidopsis (in bold). The numbers near the tree branches represent bootstrap values. Figure S7: Phylogenetic tree of PRMT proteins of Taraxacum kok-saghyz and Arabidopsis (in bold). The numbers near the tree branches represent bootstrap values. Table S1: List of T. kok-saghyz HMGs. Table S2: Primer list. Table S3: Pairs of putative duplicated HM genes in the *Tks* genome. Table S4: Tukey's pairwise test

**Author Contributions:** Conceptualization, F.P., L.L. and G.P.; methodology, F.P., L.L., A.P., E.F. and L.D.; validation, C.F. and P.F.; writing—review and editing, F.P., C.F., L.L., A.P., P.F., E.F., L.D. and G.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was partly funded by University of Sassari (FAR 2020 to AP.).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** RNAseq reads archived as part of the PRJCA000437 bioproject were downloaded from https://ngdc.cncb.ac.cn/search/?dbId=gwh&q=%20PRJCA000437 (accessed at 15 May 2022) which is included in the CNCB/NGDC genome warehouse.

**Acknowledgments:** We would like to thank USDA-ARS (Regional Plant Introduction Station, Washington State University, USA for donating the *Tks* population used in this study and National Genomics Data Center (NGDC), part of the China National Center for Bioinformation (CNCB) for making publicly available the RNA reads from Taraxacum RNAseq. Graphical abstract was created with Biorender.com.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

### *Article* β**-Cyclocitral Does Not Contribute to Singlet Oxygen-Signalling in Algae, but May Down-Regulate Chlorophyll Synthesis**

**Thomas Roach \*, Theresa Baur and Ilse Kranner**

Department of Botany, University of Innsbruck, Sternwartestraße 15, 6020 Innsbruck, Austria **\*** Correspondence<EMAIL_ADDRESS>Tel.: +43-512-507-51030

**Abstract:** Light stress signalling in algae and plants is partially orchestrated by singlet oxygen (1O2), a reactive oxygen species (ROS) that causes significant damage within the chloroplast, such as lipid peroxidation. In the vicinity of the photosystem II reaction centre, a major source of 1O2, are two β-carotene molecules that quench 1O2 to ground-state oxygen. 1O2 can oxidise β-carotene to release β-cyclocitral, which has emerged as a 1O2-mediated stress signal in the plant *Arabidopsis thaliana*. We investigated if β-cyclocitral can have similar retrograde signalling properties in the unicellular alga *Chlamydomonas reinhardtii*. Using RNA-Seq, we show that genes up-regulated in response to exogenous β-cyclocitral included *CAROTENOID CLEAVAGE DIOXYGENASE 8* (*CCD8*), while down-regulated genes included those associated with porphyrin and chlorophyll anabolism, such as tetrapyrrole-binding protein (*GUN4*), magnesium chelatases (*CHLI1*, *CHLI2*, *CHLD*, *CHLH1*), light-dependent protochlorophyllide reductase (*POR1*), copper target 1 protein (*CTH1*), and coproporphyrinogen III oxidase (*CPX1*). Down-regulation of this pathway has also been shown in β-cyclocitral-treated *A. thaliana*, indicating conservation of this signalling mechanism in plants. However, in contrast to *A. thaliana*, a very limited overlap in differential gene expression was found in β-cyclocitral-treated and 1O2-treated *C. reinhardtii*. Furthermore, exogenous treatment with βcyclocitral did not induce tolerance to 1O2. We conclude that while β-cyclocitral may down-regulate chlorophyll synthesis, it does not seem to contribute to 1O2-mediated high light stress signalling in algae.

**Keywords:** high light stress; singlet oxygen; signalling; GPX5; beta cyclocitral; acrolein; glutathione peroxidase; carbonyl; transcription

#### **1. Introduction**

Photosynthetic organisms often encounter suboptimal conditions, leading to the absorption of excess light energy and 'high light' (HL) stress. Therefore, light harvesting must be regulated, requiring acclimation to the current environment [1]. As part of acclimation, signals from the chloroplast can alter transcription in the nucleus in so-called retrograde signalling [2,3]. For example, levels of tetrapyrrole intermediates (chlorophyll precursors) in the chloroplast provide feedback cues to the nucleus during chlorophyll synthesis [4,5]. Stress signalling is also partially orchestrated by reactive oxygen species (ROS), including singlet oxygen (1O2), which is produced by energy transfer from excited chlorophyll to molecular oxygen in photosystem II (PSII) [6–8]. Singlet oxygen oxidises almost anything in its path and the 1O2 signal leaves the chloroplast in the form of down-stream reaction products. β-cyclocitral, a 1O2-derived breakdown product of β-carotene, has emerged as an aldehyde electrophile involved in 1O2 retrograde signalling of *Arabidopsis thaliana* [9]. However, the contribution of β-cyclocitral to 1O2 signalling in algae is unknown.

Other 1O2-derived molecules with potent signalling activity include α,β-unsaturated carbonyl derivatives, known as reactive carbonyl/electrophile species (RES). These are produced as a consequence of lipid peroxidation [10,11]. Thylakoid membranes are particularly

**Citation:** Roach, T.; Baur, T.; Kranner, I. β-Cyclocitral Does Not Contribute to Singlet Oxygen-Signalling in Algae, but May Down-Regulate Chlorophyll Synthesis. *Plants* **2022**, *11*, 2155. https://doi.org/10.3390/ plants11162155

Academic Editors: Milan S. Stankovic, Paula Baptista and Petronia Carillo

Received: 2 June 2022 Accepted: 15 August 2022 Published: 19 August 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

enriched in 1-linolenoyl/2-linolenoyl (di-Cl8:3) found in monogalactosyldiacylglycerol, a polar lipid that improves fluidity for membrane functionality, but is highly prone to peroxidation by 1O2. One of the most abundant RES produced by lipid peroxidation in chloroplasts due to HL stress is acrolein [10,11], which can activate a significant proportion of transcriptional changes that occur in response to 1O2 [12,13]. In *Chlamydomonas reinhardtii*, an electrophile response element (ERE)-containing bZIP transcription factor called SOR1 participates in 1O2 signalling [12,13]. SOR1 up-regulates transcription of a large suite of genes, including glutathione transferases (e.g., *GSTS1*) and an isoflavone reductase-like protein (*IRL1*), which contribute to RES-associated antioxidant defences [12,14,15].

Although they are potent signalling molecules, RES are also toxic to cells by forming Michael addition adducts with nucleophilic thiolate anions [16], such as redox-active cysteine residues of proteins. As a thiol, the antioxidant glutathione (GSH), can detoxify RES. *Chlamydomonas reinhardtii* rapidly responds to low concentrations of the RES acrolein (≤600 ppm) by increasing GSH contents, but is critically GSH-depleted at higher acrolein concentrations [13]. Another component of RES detoxification putatively includes glutathione peroxidases (GPX), whereby *GPX5* (also known as GPXh) transcription is strongly up-regulated by sub-lethal levels of 1O2 [17,18] and acrolein [13]. Moreover, overexpression of GPX5 in *C. reinhardtii* can increase tolerance to 1O2 [18].

In summary, low levels of ROS/RES can activate signalling pathways implicated in acclimation, whereas an excess ROS/RES load leads to intolerable stress. Thus, stress responses can be distinguished into eustress and distress: Eustress leads to increased stress tolerance via acclimation in which 1O2 signalling can be involved, and distress leads to loss of viability due to excess stress (e.g., high RES load) beyond a level that can be compensated for by acclimation [18–20].

Here, we explored the potential role for β-cyclocitral in 1O2-mediated signalling and inducing 1O2 tolerance in the unicellular model green alga, *C. reinhardtii*. First, we confirmed that β-cyclocitral could enter cells by measuring increased concentrations of the molecule in treated cells, and observing the concentration-dependent effect on chlorophyll fluorescence. Then, using RNA-Seq analysis, we analysed the transcriptional response of cells to β-cyclocitral and compared this to previously published data of differential gene expression induced by the photosensitizer rose bengal (RB), to reveal if elements of 1O2 signalling were activated by β-cyclocitral. Further experiments investigated if βcyclocitral can induce 1O2 tolerance, and comparisons are drawn with responses to acrolein, a 1O2-derived signal resulting from lipid peroxidation of the thylakoid membrane.

#### **2. Results**

As an aldehyde, β-cyclocitral could be cytotoxic as well as a signalling molecule. Therefore, to assess toxicity, the impact of various concentrations of β-cyclocitral on photosynthesis was probed via chlorophyll fluorescence. The maximum quantum yield of PSII (*F*v/*F*m), which is an often-used health marker of photosynthetic organisms, decreased in both wild types (WT), cell-wall containing WT-4A and cell-wall-less cw15, in response to β-cyclocitral treatment with ≥10 μL/Petri dish (Figure 1A), corresponding to ≥50,000 ppm atmospheric concentration (see methods for calculation). In contrast, NPQ that is also measured via chlorophyll fluorescence was affected at much lower β-cyclocitral concentrations (Figure 1B), including at 0.12 μL/Petri dish (Figure S1), which corresponds to the 600 ppm treatment used for the RNA-Seq analysis. Cellular concentrations of β-cyclocitral before treatment were 0.2 nmol g−<sup>1</sup> fresh weight, and increased 400 fold 2 h after exogenous treatment at 600 ppm (Figure S2). There was no difference in the influence of β-cyclocitral on NPQ or *F*v/*F*m between WT-4A and cw15 strains (Figure 1). Furthermore, the reduction of NPQ was also found in the *npq4*, *stt7,* and *npq4stt7* mutants (Figure S1), and thus was independent of LHCSR3- and STT7-mediated NPQ, which are the major NPQ mechanisms in *C. reinhardtii* [1,21]. Therefore, our results are indicative of a direct physical effect of β-cyclocitral on NPQ that occurs at the level of the thylakoid membrane, similar to an uncoupler. To test this, the activity of the pH-dependent violaxanthin cycle was measured

in response to HL. Less zeaxanthin accumulated in β-cyclocitral-treated cells (Figure S1), indicating that the proton gradient was dissipated by β-cyclocitral.

**Figure 1.** Effect of exogenous β-cyclocitral on *F*v/*F*<sup>m</sup> and NPQ in *Chlamydomonas reinhardtii.* (**A**) Measurements of *F*v/*F*m and (**B**) NPQ in 4 h HL-acclimated cultures were made after 4 h of treatment under very LL (see methods), of cell-wall-containing wild type (WT-4A; black bars) and cell-wall-less (cw15; grey bars) cultures, *n* = 3 ± SD, with distinct letters indicating significant differences (*p* < 0.05).

Aldehydes can be detoxified by GSH and associated enzymes, such as GSTS1. Thus, the response of β-cyclocitral on GSH concentrations was also measured up to 1500 ppm (0.3 μL/Petri dish), and after 4 h, no effect was observed (Figure S3).

Overall, the lack of decrease in *F*v/*F*<sup>m</sup> and lack of change in GSH concentrations indicated that β-cyclocitral was not a cytotoxic aldehyde at treatment with up to 1500 ppm, a concentration at which the RES acrolein is lethal [13]. Moreover, the significant impact of β-cyclocitral on NPQ in WT-4A and cw15 showed that this aldehyde could enter the chloroplast of *C. reinhardtii* and that the cell wall was not a hindrance to influx. Therefore, we felt confident to be able to assess the signalling properties of β-cyclocitral in exogenouslytreated cells.

One of the most 1O2-responsive genes in *C. reinhardtii* is *GPX5*, and GPXs likely have a role in mitigating HL stress by detoxifying aldehydes/RES that are produced as a consequence of 1O2 production [18]. We screened three mutants over-expressing *GPX5* in the WT-4A background [18] and found clearly elevated protein levels of GPX5 in *GPXHOX-11* and *GPX5OX-14*, relative to WT-4A, in low light (LL)-treated cells (Figure 2A). This led to the selection of *GPXHOX-11* for further experiments to see how β-cyclocitral, acrolein, and HL influence GPX5 levels, and how elevated levels of GPX5 affect the response to these treatments and subsequent tolerance to 1O2. Levels of GPX5 increased in WT-4A cells after 4 h of HL, indicating that this treatment induced 1O2 production. However, 4 h treatments with 600 ppm of β-cyclocitral or 600 ppm of acrolein did not increase GPX5 levels (Figure 2B). After 4 h of HL, WT cells increased the GSH contents, indicating that cells were under mild oxidative stress, whereas in *GPXHOX-11,* no change in the GSH contents occurred, indicating that cells were under less stress. In support of this, after 4 h of HL, *GPXHOX-11* accumulated less aldehydes (propanal and hexanal) and RES (acrolein and 4-hydroxynonenal) than WT-4A, whereas the fold change of β-cyclocitral was equally low in both genotypes (Figure 3). After a 4 h treatment with 600 ppm of acrolein, both WT and *GPXHOX-11* accumulated GSH (Figure S3).

**Figure 2.** Protein levels of GPX5 in WT-4A and *GPX5*-overexpressor (*GPXHOX)* mutants under LL and in response to HL, β-cyclocitral and acrolein. (**A**) Three *GPXHOX* lines and WT-4A were analysed for GPX5 protein levels under LL. Shown above is the Ponceau-stained membrane for loading control. (**B**) The effect of LL, 4 h with HL, and 4 h with 600 ppm of β-cyclocitral (β-cc) or with 600 ppm of acrolein (4 h Acro.) under LL, on GPX5 levels. The D1 reaction centre of photosystem II (PsbA) was used for loading control.

**Figure 3.** Influence of GPX5 overexpression on aldehyde/RES accumulation. Change in RES concentrations in response to HL (600 μmol photons m−<sup>2</sup> s<sup>−</sup>1) for 4 h, relative to levels before treatment at LL (50 μmol photons m−<sup>2</sup> s<sup>−</sup>1), in wild-type (black) and *GPXHOX-11* (white). MDA: malondialdehyde, βcc: β-cyclocitral, 4HNE: 4-hydroxynonenal. *n* = 4 ± SD with \* denoting significant differences between genotypes (*p* < 0.05).

A function of 1O2-related signalling is partly towards increasing tolerance of this ROS [18]. Therefore, the influence of exogenous β-cyclocitral and acrolein treatments under LL (as potential components of 1O2-mediated signalling) on the tolerance to 1O2 was measured. The influence of HL stress, which increases endogenous 1O2 levels, was also included as pre-treatment before testing 1O2 tolerance, and all comparisons were made to LL-treated 'control' cells. Tolerance to 1O2 was tested by incubating cells with the photosensitizer RB. Since the amount of 1O2 that RB produces is dependent on the degree of photoexcitation (i.e., intensity and duration of light treatment) and RB concentration, various treatments were conducted to test cell tolerance. These constituted 4 μM, 7 μM, and 10 μM of RB, either for 10 min at 250 μmol photons m−<sup>2</sup> s−<sup>1</sup> to provide a short 1O2 shock, or 24 h at 50 μmol photons m−<sup>2</sup> s−<sup>1</sup> to test longer term endurance. The 24 h endurance test was more severe and led to less cell survival of control cells, under which *GPXHOX-11* was significantly less affected than WT-4A (Figure 4), with *p* < 0.05 when comparing genotype as a factor with MANOVA across all RB concentrations. Pre-treatment with HL increased the tolerance of WT-4A to the 1O2 shock treatment with 7 μM and 10 μM of RB (Figure 4A), and the 24 h treatment with 4 μM (Figure 4B), while affecting *GPXHOX-11* less. This resulted in HL-treated WT-4A and *GPXHOX-11* having a similar 1O2 tolerance (Figure 4). Relative to control cells, pre-treatment with 600 ppm of β-cyclocitral had no impact on 1O2 tolerance in either genotype under all treatments (average *p* = 0.895), whereas pre-treatment with 600 ppm of acrolein increased the tolerance of both genotypes (average *p* = 0.034). In summary, acrolein and HL induced tolerance to severe 1O2 stress, whereas β-cyclocitral did not significantly affect tolerance.

**Figure 4.** Influence of *GPX5*-overexpression and effect of HL, β-cyclocitral, and acrolein pretreatments on tolerance to 1O2. Tolerance to 1O2 was tested by exposing WT 4A or *GPXOX-11* cultures to (**A**) 250 μmol photons m−<sup>2</sup> s−<sup>1</sup> for 10 min or (**B**) 50 μmol photons m−<sup>2</sup> s−<sup>1</sup> for 24 h, in a medium that initially contained 4 μM, 7 μM, or 10 μM of rose bengal (RB). Before treatment with RB, cells were from LL (control; white) or pre-treated with 4 h of HL (black), 600 ppm of acrolein (blue) or 600 ppm of β-cyclocitral (red) under LL. To reveal tolerance to 1O2, cell density was monitored for 24 h after RB treatment by turbidity of culture at 650 nm and is shown relative to turbidity of the respective pre-treated cultures not treated with RB. Bars labelled with different letters within the same RB treatment denote significant differences from each other (Tukey's post hoc after arcsine transformation, *p* < 0.05, *n* = 3 ± SD).

To investigate if β-cyclocitral could contribute to 1O2-mediated signalling in *C. reinhardtii*, an RNA-Seq analysis of cells treated for 2 h with 600 ppm of β-cyclocitral under LL (Table S1) was conducted and compared to the transcriptional response of cells treated with the 1O2-producing photosensitizer RB at 1 μM under LL (data from [13]). In response to the β-cyclocitral treatment, only six genes were significantly up-regulated, and 57 genes were down-regulated, when considering a fold change of >2 and modified *t*-test *p* values of <0.01 (Table S1). Of the genes down-regulated by β-cyclocitral, 18 and 6 were significantly downregulated and up-regulated, respectively, by RB (Table S1; Figure S4). Differentially expressed genes associated with carotenoid metabolism include carotenoid cleavage dioxygenase 8 (*Cre08.g365851*), up-regulated five-fold, and β-carotene ketolase (*BKT1*), up-regulated two-fold, but not significantly (Table S1). Of all significantly down-regulated genes, the only ontological group with >1 hit was 'porphyrin and chlorophyll metabolism' with 10 hits: *Cre01.g015350* (*POR1*), *Cre01.g050950*, *Cre02.g085450* (*CPX1*), *Cre05.g242000* (*CHLD*), *Cre05.g246800* (*GUN4*), Cre06*.g306300* (*CHLI1*), *Cre07.g325500* (*CHLH1*), *Cre09.g396300* (*PPX1*), *Cre12.g510050* (*CTH1*), *Cre12.g510800* (*CHLI2*), and *Cre16.g663900*. Collectively, these genes covered several steps of chlorophyll anabolism (Figure S5), but the overall overlap of differential gene expression with RB-treated cells was low, as shown by a R2 linear correlation of <0.01 when including all genes, which increases to 0.15, considering only the 63 genes with significantly altered expression (Figure 5). In comparison, this contrasts with the much tighter relationship between differential gene expression shared between acrolein-treated and RB-treated cells [13], which has an R2 linear correlation of 0.34 and 0.70 when considering expression of total genes and only significantly affected genes, respectively (Figure 5).

**Figure 5.** Correlation of differential gene expression induced by treatment with 1O2 and acrolein (blue) or 1O2 and β-cyclocitral (red). Levels of mRNA are expressed as log2 fold changes (log2FC), relative to mock-treated cells, calculated from RNA-Seq analyses; *n* = 3 for treatments and controls. Treatments with acrolein and β-cyclocitral were at 600 ppm and 1O2 was induced by rose bengal (RB; 1 μM under growth light at 50 μmol quanta m−<sup>2</sup> s<sup>−</sup>1). Continuous and dashed lines of best fit consider all and only significantly affected genes, respectively. See Table S1 for expression levels of each gene. RNA-Seq data for acrolein and RB treatments adapted with permission from [13] 2018, Copyright Elsevier.

#### **3. Discussion**

Ten years ago, the discovery that β-cyclocitral in *A. thaliana* effects transcription of several genes know to be affected by 1O2 [9] made a coherent link between HL-induced 1O2 production and ROS-associated retrograde signalling. Subsequently, it was found that β-cyclocitral functions up-stream of MBS1 [22], a zinc finger protein that regulates 1O2-dependent gene expression, not only in *A. thaliana* but also in *C. reinhardtii* [23]. However, the fact of whether β-cyclocitral actually has a role in 1O2 signalling in alga remained unknown. Since then, other RES (i.e., acrolein), related to lipid peroxidation rather than carotenoid cleavage, emerged as retrograde signals acting in 1O2-mediated stress acclimation [13]. Here, we investigated how β-cyclocitral modulates the physiology and transcription in *C. reinhardtii* and made comparisons with transcriptional responses of cells to the photosensitizer RB and the RES acrolein.

The very minor effect of β-cyclocitral on *F*v/*F*<sup>m</sup> showed how tolerant cells were of this molecule. For example, decreases in *F*v/*F*<sup>m</sup> occurred at >1000 fold concentrations compared to the effects of acrolein (Figure 1A; [13]). The concentration of β-cyclocitral even in very light-stressed *C. reinhardtii* has never been measured at >1 nmol g−<sup>1</sup> fresh weight [13,24], which is below the cellular concentrations after exogenous treatment with 600 ppm (Figure S2). Therefore, our data supports that, unlike acrolein, β-cyclocitral does not build up to toxic concentrations in light-stressed cells. Despite the insensitivity of *F*v/*F*<sup>m</sup> to β-cyclocitral, NPQ was affected at much lower concentrations (Figures 1B and S1), confirming that β-cyclocitral was able to enter cells, which otherwise may have contributed to tolerance. Inhibition of NPQ in various NPQ mutants indicates that β-cyclocitral directly affected an over-riding NPQ mechanism, such as the requirement of low luminal pH. This was indeed shown by the lower accumulation of zeaxanthin under HL, which is a process requiring a low luminal pH for violaxanthin de-epoxidase activity, in β-cyclocitral-treated cells (Figure S1). Therefore, we suggest that β-cyclocitral may act as an uncoupler of the thylakoid membrane potential, comparable to the inhibitory activity of structurally similar monoterpene ketones, such as pulegone, on respiration [25]. The relevance of this observation is that NPQ protects from 1O2 production under HL [1,24], and thus exogenous β-cyclocitral may increase 1O2 production and may confound observations of β-cyclocitral involvement in 1O2 signalling under HL. Here, exogenous treatments for the RNA-Seq were conducted under very LL (2 μmol photons m−<sup>2</sup> s<sup>−</sup>1), thus unaffected by lowered NPQ.

Defence against 1O2 requires many enzymes that mitigate lipid peroxidation, including GPX5 [18]. In mammalian cells, GSH is the typical GPX substrate to break down H2O2, whereas in *C. reinhardtii,* a thioredoxin is the reductant of GPX5 that has a close association with 1O2 stress [26]. In WT cells, HL induced accumulation of GPX5 (Figure 2B), alongside higher levels of RES, which were attenuated in *GPXHOX-11* (Figure 3), a *GPX5* overexpressing mutant with elevated GPX5 levels (Figure 2). This supports that GPX5 has a role in metabolising HL-induced RES production, and explains why HL stress is associated with elevated GPX5 levels [26,27]. *GPXHOX-11* also possessed significantly elevated tolerance to long-term 1O2 treatment with the photosensitizer RB (Figure 4B), in agreement with results from Ledford et al. [18]. Comparing the influence of pre-treatments on 1O2 tolerance, acrolein was able to induce tolerance of WT and *GPXHOX-11*, whereas β-cyclocitral could not (Figure 4). B-cyclocitral did not induce GSH synthesis, whereas acrolein and HL did (Figure S3), as also previously shown [13]. Thus, enhanced 1O2 tolerance can be partially attributed to elevated GSH and GPX5 levels.

A molecule involved in stress signalling would be expected to increase in concentration in response to the relevant stress. The amounts of β-cyclocitral in *C. reinhardtii* were ≤1 nmol g−<sup>1</sup> fresh weight (Figure S2), and not affected by HL stress (Figures S2 and 3), also in agreement with previous data [13,24]. In comparison, levels of RES increased in HL-stressed cells, similar to after treatment with RB [13], supporting that 1O2 production under HL is involved in RES production, but hardly with β-cyclocitral production. Acrolein has received attention in the field of redox biology for its high electrophilic nature and high endogenous levels of >5 nmol g−<sup>1</sup> fresh weight in stressed plants and algae alike [13,28,29]. Previously, we showed that exogenous acrolein treatments induce a 'eustress' (i.e., acclimation) response by up-regulating thiol-disulfide-dependent defence mechanisms required for tolerating 1O2 [13]. Of note, around half of global gene expression (up and down) occurring in response to 600 ppm of acrolein, the dose that induced a eustress response with highest tolerance to 1O2, was shared with the gene regulation in response to RB [13]. While we did not find a similar transcriptional response to β-cyclocitral (Figure 5), in line with a lack of inducing 1O2 tolerance (Figure 4), there was evidence that β-cyclocitral may have some signalling properties in C*. reinhardtii*. Collectively down-regulated genes (Table S1) covered many steps of chlorophyll anabolism (Figure S5), such as porphobilinogen deaminase/HemC (*Cre16.g663900.t1.2*), coproporphyrinogen III oxidase (*CPX1*), and protoporphyrinogen oxidase (*PPX1*), which are involved in early steps of porphyrin synthesis, as well as genes coding for proteins that insert Mg2+ into protoporphyrin, including tetrapyrrole-binding protein (*GUN4*) and magnesium chelatase (*CHLI1*, *CHLI2*, *CHLD*, *CHLH1*) to form Mg-protoporphyrin IX (MgP), the first dedicated intermediate of the