- 1Department of Botany and Plant Sciences, University of California, Riverside, Riverside, CA, United States
- 2U.S. Horticultural Research Laboratory, Agricultural Research Service, United States Department of Agriculture, Ft. Pierce, FL, United States
Interspecific hybridization is a common and effective strategy for producing disease resilient citrus cultivars, including those with tolerance to Huanglongbing (HLB) disease. Several HLB-tolerant cultivars have been developed through hybridization of mandarins (Citrus reticulata) with their wild relative Poncirus trifoliata. One such cultivar, ‘US-897’, exhibits robust tolerance to the bacteria causing HLB disease, Candidatus Liberibacter asiaticus (CLas). To explore the genetic architecture of the early transcriptional response to Candidatus Liberibacter asiaticus (CLas) infection in ‘US-897’, we performed transcriptomic analysis of the hybrid and its parents, ‘Cleopatra’ (C. reticulata) and ‘Flying Dragon’ (P. trifoliata). A haplotype-resolved genome for ‘US-897’ was generated using PacBio HiFi sequencing reads to support quantification of the expression of both the Citrus and Poncirus alleles. By profiling gene expression in this parent-offspring trio, we were able to determine the mode of inheritance for genes differentially expressed between parents (‘Cleopatra’ and ‘Flying Dragon’) and their interspecific hybrid (‘US-897’), with the majority genes exhibiting non-additive patterns of gene expression inheritance. Additionally, analysis of allele-specific expression in the hybrid ‘US-897’ revealed the contribution of cis- versus trans-acting regulatory variants on genes with additive and non-additive modes of inheritance. A strong correlation between differential expression between parents and allele-specific expression in ‘US-897’ suggests that cis-regulatory variation is a significant source of expression divergence between species. Finally, genes responsive to infection with CLas were identified to explore how gene regulation associated with tolerance to HLB was rewired between Citrus and its relative Poncirus.
Introduction
Interspecific hybridization can improve disease resilience by introducing beneficial alleles from the relatives of cultivated crop species. In long-lived perennial tree crops, like those in the genus Citrus, interspecific hybridization has led to the development of many commercial cultivars with enhanced disease resistance and improved productivity in disease endemic regions. A current major challenge to citrus production is the bacterial disease Huanglongbing (HLB). The most common cause of HLB in the United States is the phloem-restricted bacterial species Candidatus Liberibacter asiaticus (CLas) which infects all commercial citrus types and has major impacts on fruit production (Bové, 2006). CLas is vectored by the Asian citrus psyllid (Diaphorina citri), which first arrived in Florida in 1998 (Hall, 2008). HLB has devastated the citrus industry in Florida (Court and Rahmani, 2017) and genetic solutions, including the development of disease resilient cultivars, are one of the many tools being used to combat HLB disease (Albrecht and Bowman, 2012a; Bowman and McCollum, 2015; Albrecht et al., 2020; Bowman and Albrecht, 2020; Peng et al., 2020).
Poncirus trifoliata is a relative of species in the genus Citrus with tolerance to HLB (Albrecht and Bowman, 2012a). Several Citrus x Poncirus hybrids have been shown to have levels of tolerance to CLas similar to P. trifoliata (Albrecht and Bowman, 2011, 2012a; Bowman and McCollum, 2015), including reduced bacterial titers and limited development of disease symptoms. One hybrid, ‘US-897’, developed from a cross between ‘Cleopatra’ mandarin (C. reticulata) and ‘Flying Dragon’, (P. trifoliata) exhibits robust tolerance to CLas infection as an ungrafted tree (Albrecht and Bowman, 2011). Although the hybrid ‘US-897’ can still be infected by the bacterium, it only exhibits subtle symptoms of infection in both greenhouse studies and field evaluations (Albrecht and Bowman, 2011). In contrast, CLas titer rapidly increases in ‘Cleopatra’, the susceptible maternal parent of ‘US-897’, after infection and associated disease symptoms include chlorosis, blotchy mottling of leaves, and reduced overall size. The strong resilience of ‘US-897’ to CLas provides an opportunity to characterize the basis of HLB tolerance in an interspecific hybrid.
Improved hybrid performance can result from genetic interactions between parental genomes, including those that may alter patterns of gene expression (Chen, 2013). In the analysis of parent-offspring trios, gene expression can be classified into additive effects, where hybrid expression is intermediate between the two parents, and non-additive effects such as expression-level dominance and transgressive expression including over- and under-dominance (Rapp et al., 2009). There has been a long debate over the genetic architecture underlying hybrid traits (Schwartz and Laughner, 1969; Semel et al., 2006; Krieger et al., 2010; Seymour et al., 2016), and the relative contribution of loci with dominant (Hua et al., 2003; Garcia et al., 2008; Xiao et al., 2021) or overdominant (Rédei, 1962; Vrebalov et al., 2002; Krieger et al., 2010; Guo et al., 2014) contributions to hybrid vigor (Paril et al., 2024). Examples of loci with dominant and overdominant effects on hybrid vigor have been identified across many plant species (Semel et al., 2006; Li et al., 2008; Hashimoto et al., 2021). Overall, the contribution of dominant and overdominant loci to hybrid vigor seems to vary by species and mating system (Huang et al., 2016; Tian et al., 2019; Hashimoto et al., 2021; Xiao et al., 2021). Ultimately, the genetic architecture and patterns of inheritance are important for understanding the regulatory mechanisms underlying hybrid phenotypes, including gene expression.
The genetic architecture of gene regulatory variation in hybrids includes whether genes are cis- or trans-regulated and can provide insight into the mechanistic basis of hybrid phenotypes (Rockman and Kruglyak, 2006). Cis-regulatory variants alter expression of a gene by disrupting regulatory sequences in close proximity to that gene and may include genetic variation in promoters and enhancers. In contrast, trans-regulatory variants perturb the expression of a gene by altering the activity of regulatory molecules, such as transcription factors and regulatory RNAs (Hill et al., 2021). Unlike cis-regulatory variants, trans-regulatory variants may not be located in close physical proximity to their target. Profiling gene expression in parent-offspring trios is commonly used to separate cis- and trans-regulatory effects on expression because both parental alleles, and any associated cis-regulatory variants, experience a shared trans-regulatory environment (Wittkopp et al., 2004; Goncalves et al., 2012; Pereira and Kuzmin, 2024). Cis-regulated genes will exhibit allele-specific expression (ASE) in the hybrid, while trans-regulatory variants will alter the expression of both alleles concordantly (Hill et al., 2021). Dissecting the genetic architecture of gene expression can connect expression divergence to phenotypic divergence and can provide valuable insight into the mechanisms underlying hybrid phenotypes.
Here, we assembled a high-quality, haplotype-resolved reference genome for the HLB-tolerant hybrid ‘US-897’ and performed transcriptomic analyses of the hybrid and its parents, ‘Cleopatra’ (C. reticulata) and ‘Flying Dragon’ (P. trifoliata), to investigate the genetic architecture of gene expression regulation, including in response to CLas infection. Haplotype phasing of the ‘Cleopatra’ and ‘Flying Dragon’ chromosomes in the diploid genome of ‘US-897’ enabled quantification of allele-specific expression (ASE) through accurate assignment of RNA-sequencing reads to either the Citrus or Poncirus haplotype, which is superior to relying on a haploid reference genome (Wu et al., 2023). By profiling gene expression in this parent-offspring trio, we were able to determine the mode of inheritance for genes differentially expressed between parents (‘Cleopatra’ and ‘Flying Dragon’) and their interspecific hybrid (‘US-897’). Additionally, analysis of allele-specific expression in the hybrid ‘US-897’ revealed the contribution of cis- versus trans-acting regulatory variants on genes with additive, dominant, and transgressive inheritance. Finally, genes responsive to infection with CLas were identified to explore the regulatory architecture of pathogen response.
Results
Assembly of a diploid reference genome of F1 hybrid ‘US-897’
To quantify the expression of both the Citrus and Poncirus alleles in the interspecific hybrid, a haplotype-resolved chromosome-level genome of ‘US-897’ was assembled from 65.57 Gb of PacBio HiFi sequencing reads (108.75X coverage; mean read length = 15.68 Kb). Parent-specific k-mers (k=19), or “hap-mers”, were identified from whole-genome sequence data for ‘Cleopatra’ mandarin (C. reticulata) and ‘Flying Dragon’ (P. trifoliata) and used to phase the Citrus and Poncirus haplotypes in ‘US-897’ (Wu et al., 2018; Peng et al., 2020; Rhie et al., 2020; Rautiainen et al., 2023). For both haplotypes, the nine largest scaffolds had BUSCO scores of 99.0%, and captured 99.2% and 99.8% of parental hap-mers, an indication of assembly completeness (Supplementary Table S1). There was no evidence of haplotype switching within scaffolds, confirming that each haplotype was correctly phased in the diploid assembly of ‘US-897’ (Figure 1A). Haplotype 1 scaffolds were inherited from ‘Flying Dragon’ (P. trifoliata) (282.4 Mb, N50 = 29.9 Mb) and haplotype 2 scaffolds from ‘Cleopatra’ (C. reticulata) (321.4 Mb, N50 = 34.1 Mb) (Figure 1A; Supplementary Table S1). A similar number of genes were annotated in each haplotype (Haplotype 1: 33,431; Haplotype 2: 31,602 genes) and, based on BUSCO scores, the annotations included a majority of single copy orthologous genes used to benchmark genome annotations (Haplotype 1: 98.8%; Haplotype 2: 98.4%) (Supplementary Table S1). Orthology between the Citrus and Poncirus gene models was inferred based on protein sequence identity and collinearity using GeneTribe (Chen et al., 2020). This identified 22,419 reciprocal-best hits that were used for comparisons of allele-specific expression, along with 13,439 orthologous genes (single-best hits) and 6,756 singletons unique to one of the two parental haplotypes (Figure 1B). Analysis of allele-specific expression requires identification of gene orthology and this set of 22,419 genes were used to dissect the regulatory architecture of gene expression divergence between Citrus and Poncirus.

Figure 1. A haplotype-resolved assembly of the interspecific hybrid ‘US-897’. (A) The abundance of parent-specific “hap-mers” (k=19) in each haplotype of the ‘US-897’ diploid reference genome. Hap-mers were identified from ‘Cleopatra’ and ‘Flying Dragon’ whole-genome sequencing reads (Wu et al., 2018; Peng et al., 2020). (B) Classification of genes based on their sequence homology and collinearity between Citrus and Poncirus haplotypes of ‘US-897’. Reciprocal best hits (RBH) are allelic between haplotypes, single best hits (SBH) are orthologous sequences that are not reciprocal best hits, and singletons are genes that are unique to a haplotype.
Expression level dominance is pervasive in ‘US-897’
To better understand tolerance to HLB disease in ‘US-897’, we profiled the early transcriptional response of ‘Cleopatra’, ‘Flying Dragon’, and ‘US-897’ to CLas infection. The three genotypes were each infested by Asian citrus psyllids reared either on CLas-infected or clean citron plants (C. medica). Infection with CLas was confirmed for five ‘Cleopatra’, two ‘Flying Dragon’, and five ‘US-897’ individuals (Supplementary Table S2). To capture an early transcriptional response to CLas infection, leaf samples were collected from newly emerged branches 18 days after bud initiation. Previous work has shown early cellular and biochemical responses to CLas infection, including reactive oxygen species (ROS) production and callose deposition, occur as early as 15 days after inoculation (Ma et al., 2022). It is important to note that this sampling scheme captures initial host responses and the relationship between initial response and the robust tolerance of ‘US-897’ to Huanglongbing disease is unknown. An average of 28.3 million paired-end RNA-sequencing reads were generated per sample (Supplementary Table S3). RNA-seq reads were aligned to the haplotype-resolved genome of ‘US-897’ and the level of expression of the Citrus and Poncirus alleles of each gene were quantified using Salmon (Patro et al., 2017). In addition to allele-specific expression, levels of overall gene expression were determined by summing read counts for the Citrus and Poncirus alleles of each orthologous gene. Genes with low levels of expression were removed (CPM less than 2) and the remaining set of 20,981 genes were the basis of differential gene expression analysis.
The analysis of gene expression in parent-offspring trios is essential for determining patterns of gene expression inheritance. Whether a gene exhibits additive, dominant, or transgressive patterns of inheritance in an F1 hybrid compared to its two parents can provide mechanistic insights into gene expression divergence between parental alleles. A linear model was constructed to explain overall gene expression by a genotype effect (‘Cleopatra’, ‘Flying Dragon’, ‘US-897’, and mid-parent value), a treatment effect (CLas-inoculated, mock-inoculated), and the interaction between genotype and treatment. A total of 9,522 differentially expressed genes (DEGs) with a significant genotype effect were identified (‘Cleopatra’ versus ‘Flying Dragon’, ‘Cleopatra’ or ‘Flying Dragon’ versus ‘US-897’, and ‘US-897’ versus estimated mid-parent value) (Figure 2A, FDR adjusted p-value < 0.05). A majority of these genes (8,807) were differentially expressed between the parental genotypes (Figure 2A). For genes with additive expression the level of expression in the hybrid is similar to the average of parental allelic expression, or the mid-parent expression value. Comparison of gene expression between ‘US-897’ and the mid-parent value revealed that gene expression patterns in the hybrid deviated from additivity for a set of 784 genes (Figure 2A).

Figure 2. Classification of patterns of gene expression inheritance in ‘Cleopatra’, ‘US-897’, and ‘Flying Dragon’. (A) Diagram of the number of differentially expressed genes (DEGs) with a significant genotype effect in comparisons of: ‘Cleopatra’ versus ‘Flying Dragon’, ‘Cleopatra’ or ‘Flying Dragon’ versus ‘US-897’, and ‘US-897’ versus estimated mid-parent value. Solid black lines indicate comparison, with the total number of DEGs per comparison noted in the center of each line and the number of upregulated DEGs per comparison noted adjacent to each genotype. Percentages are based on the total number of tested genes (n = 20,981). (B) Categorization of patterns of gene expression inheritance. The focal dot-plot compares the log2 fold-change in gene expression between ‘US-897’ and its parents, with points colored by the subcategory of expression inheritance. The expression pattern of each category is shown in the 12 inset plots. (C) The 12 subcategories of gene expression inheritance were grouped into five major classifications: transgressive up in ‘US-897’, transgressive down in ‘US-897’, additive, dominant for ‘Flying Dragon’ allele, and dominant for ‘Cleopatra’ allele.
The relative abundance of gene expression across all three genotypes was used to further categorize DEGs with significant genotype effects based on inheritance patterns (Rapp et al., 2009; Almeida-Silva et al., 2023). Genes with additive gene expression in the hybrid will be intermediate between the parental levels of expression (similar to the mid-parent value). Genes with non-additive gene expression in the hybrid will be biased towards expression of one of the two parental alleles (dominant) or is significantly higher or lower than expression of both parental alleles (transgressive). Gene expression patterns in each of these three categories can be further partitioned based on which parental allele is expressed at a higher level. DEGs were grouped into five categories (1 additive, 2 dominant, 2 transgressive) and 12 subcategories based on relative expression differences between the parental alleles (Figures 2B, C). For example, subcategories 1 and 12 both include genes with additive inheritance, but category 1 contains genes with reduced expression of the ‘Flying Dragon’ allele and category 12 includes genes with increased expression of the ‘Flying Dragon’ allele (Figure 2C). The majority of DEGs (80.83%) display expression level dominance, indicating that one of the two parental alleles is commonly favored in the hybrid ‘US-897’ (Figure 2C). Dominance of the maternal ‘Cleopatra’ allele was more common (43.45% of DEGs) than dominance of the paternal ‘Flying-Dragon’ allele (37.38% of DEGs) (Figure 2C). In Arabidopsis thaliana hybrids, expression level dominance in F1 hybrids was also common, affecting 90% of DEGs (Yuan et al., 2023). The preference for maternal expression bias has also been observed in other studies of parent-offspring trios hybrids (Videvall et al., 2016). Examples of DEGs with transgressive patterns of expression were rare (3.18%). Comparison of gene expression patterns in this parent-offspring trio revealed genes with additive, dominant, and transgressive expression patterns, but that expression level dominance was the most pervasive pattern of gene expression in ‘US-897’.
Genes with additive and dominant patterns of gene expression are regulated in cis
To better understand the evolution of gene regulation it is crucial to determine whether gene expression is controlled by cis-acting or trans-acting variants. Classification of cis- and trans-regulatory effects on gene expression can be determined by comparing the expression of parental alleles in an F1 hybrid. If differential expression between parental species (i.e. C. reticulata and P. trifoliata) is caused by a cis-regulatory variant, allele-specific expression of this gene is expected in the F1 hybrid. In contrast, a trans-regulatory variant would cause differential gene expression between parental species, but concordant expression of both parental alleles in the hybrid. Cis-regulatory variants are generally considered to be major contributors to expression divergence between different species, while trans-regulatory variants are more commonly associated with expression differences within species (Emerson et al., 2010; Metzger et al., 2017; Signor and Nuzhdin, 2018; Zhang and Emerson, 2019; Coolon et al., 2014). To determine if gene expression is controlled in cis or trans we compared DEGs with significant genotype effects to allele-specific expression (ASE) of parental alleles in ‘US-897’. Genes with ASE in ‘US-897’ were identified using the same set of 22,419 orthologous genes between Citrus and Poncirus haplotypes. After removing genes with low expression, monoallelic expression, and insufficient genetic variation (a limitation in the analysis of ASE), ASE could be evaluated for 15,549 genes. Allele-specific expression was common in ‘US-897’ (9,189 genes, FDR-corrected p < 0.05), indicating that cis-regulatory variation has a major role in regulating hybrid gene expression patterns. On average, there was a 97% difference in expression between parental alleles for genes with significant ASE (log2(Haplotype 1 (P. trifoliata)/Haplotype 2 (C. reticulata)) = 0.98). Notably, allele-specific expression in ‘US-897’ explains 41.6% of differential gene expression variation between ‘Cleopatra’ and ‘Flying Dragon’ (Figure 3A; Spearman’s rho = 0.74, R2 = 0.416). Parental haplotypes have diverged for ~9 million years (Wu et al., 2018; Huang et al., 2023b), and our results are consistent with the role of cis-regulatory variation in divergence of gene expression patterns.

Figure 3. The role of cis-regulatory variation in expression divergence between Citrus reticulata and Poncirus trifoliata. (A) The relationship between ASE (log2(FD/CLEO)Allelic) and expression divergence between parental genotypes ‘Flying Dragon’ and ‘Cleopatra’ (log2(FD/CLEO)Parental). (B) The magnitude of ASE for genes with different patterns of gene expression inheritance. Categories 3, 5, 6, and 10 are not shown because there were no or few genes with ASE.
Next, we explored if cis-regulated genes (i.e. ASE genes in ‘US-897’) corresponded to additive, dominant, or transgressive patterns of expression inheritance. On average, the allele expressed more highly in the hybrid was consistent with parental expression patterns. For example, additively expressed genes with higher expression in ‘Cleopatra’ (Category 1) exhibit allele expression bias towards the ‘Cleopatra’ allele in ‘US-897’ (Figure 3B). This trend is consistent for all patterns of additive and dominant gene expression categories (Figure 3B). Allele-specific expression was not detected for genes with transgressive expression in the F1 hybrid consistent with trans-regulation of this gene set (Figure 3B). Additionally, genes with the greatest magnitude of allele-specific expression are additively inherited, perhaps with cis-regulatory variants of large effect. The abundance of genes with additive and dominant patterns of expression inheritance and their allelic dynamics in ‘US-897’ is consistent with the major contribution of cis-regulatory variation to gene expression divergence between Citrus and Poncirus.
Genes responsive to CLas infection have reduced signature of cis-regulation
The relative importance of cis- and trans-regulation to environmental response is still unclear. There is evidence of both cis- and trans-acting regulation of gene expression in response to environmental signals (Snoek et al., 2017; Reynoso et al., 2019; Zeng et al., 2019; Zhou et al., 2021), although some suspect that trans-regulators may be more critical for environmental response (Signor and Nuzhdin, 2018). To gain insight into the role of cis- versus trans-regulation of genes conferring HLB-tolerance in ‘US-897’, we first identified genes with significant treatment effects and genotype x treatment interactions. Only 629 genes were identified with a significant treatment effect (FDR corrected p < 0.05, Supplementary Table S4). This set of genes is responsive to CLas infection regardless of genotype. Of these 629 genes, 245 were upregulated, and 384 were downregulated (Supplementary Figure S1A). The few genes responding to CLas infection here may be due to the early sampling only 18 days after bud initiation. Previous gene expression analysis comparing ‘Cleopatra’ to ‘US-897’ identified many more differentially expressed genes, but plants were sampled 7 to 8 months after infection (Albrecht and Bowman, 2012b). While early transcriptomic changes are important for understanding the initial stages of host-pathogen interactions, they do not explain the full suite of mechanisms underlying tolerance to Huanglongbing disease. Gene ontology enrichment analysis revealed that up-regulated genes were enriched for several biological processes related to regulation of the cell cycle, while down-regulated genes were enriched for processes related to cutin biosynthesis, cell wall organization, and lipid metabolism (Supplementary Table S5; Fisher’s exact test, p < 0.001). A larger group of genes with a significant genotype x treatment interaction were identified. Of the 1,071 genes with significant genotype-specific responses to HLB disease (FDR corrected p < 0.05), a majority (972) exhibited differential response to CLas infection between the two parental genotypes, ‘Cleopatra’ and ‘Flying-Dragon’. A smaller number of genes (490) were differentially responsive to CLas between one parent, ‘Flying-Dragon’, and ‘US-897’ (Figure 4A; FDR corrected p < 0.05). There was some overlap between genes with significant genotype x treatment interactions between contrasts. Surprisingly, no genes were differentially responsive to CLas between the second parent, ‘Cleopatra’ and ‘US-897’. The two parents, ‘Flying Dragon’ and ‘Cleopatra’, are responding differently to infection with CLas, at least in the early stages of infection. Few genes are differentially responsive to CLas in the hybrid compared to its two parents, suggesting that much of the parental response to infection was inherited in the hybrid.

Figure 4. Classification of patterns of gene expression inheritance in ‘Cleopatra’, ‘US-897’, and ‘Flying Dragon’ for genes with significant genotype x treatment interactions. (A) Diagram of the number of differentially expressed genes (DEGs) in comparisons of: ‘Cleopatra’ versus ‘Flying Dragon’, ‘Cleopatra’ or ‘Flying Dragon’ versus ‘US-897’, and ‘US-897’ versus estimated mid-parent value. Solid black lines indicate comparisons, with the total number of DEGs noted in the center of each line and the number of upregulated DEGs noted adjacent to each genotype. Percentages are based on the total number of test genes (n = 20,981). (B) Categorization of gene expression patterns relative to expression in ‘US-897’ versus parental genotypes for genes with significant genotype x treatment interaction (n=1,071). The dot-plot compares the log2 fold-change in gene expression between infected samples of ‘US-897’ and its two parents, with points colored by quadrant. (C) The relationship between Δ ASE in ‘US-897’ (the allelic response to CLas infection) (log2(FD/CLEO)Allelic) and expression divergence between parental genotypes ‘Flying Dragon’ and ‘Cleopatra’ (log2(FD/CLEO)Parental). Only genes with significant genotype x treatment effects that could be tested for ASE are included (n=975). (D) The response of each allele in ‘US-897’ to pathogen infection for genes in each of the quadrants in 4B. Only genes with significant genotype x treatment effects that could be tested for ASE are included (n=975).
To assess the contribution of cis-regulation to CLas responsive genes we quantified the response in allele-specific expression in ‘US-897’ (ΔASE) compared to expression divergence between the infected parents. Of the 1,071 genes with significant genotype x treatment effects, 975 could be assessed for allele-specific expression in ‘US-897’. For these genes there is a moderate correlation between ΔASE and differential expression between infected parents (Figure 4C, Spearman’s rho = 0.48, R2 = 0.264) with Δ ASE explaining 26.4% of variation in divergent expression between ‘Flying Dragon’ and ‘Cleopatra’. We also examined the relationship between ΔASE and differential expression between infected parents for 536 of the 629 genes with significant treatment effects. For these genes, ΔASE explained an even smaller proportion of variation in divergent expression between parents (R2 = 0.134; Supplementary Figure S1C). This suggests that cis-regulatory variants are also a component of divergence in disease response, but that the signature of cis-regulation is reduced.
To further investigate cis-regulation of genes responding to CLas infection, the allelic response to pathogen infection was evaluated in ‘US-897’. Genes were partitioned into quadrants that were based on whether they were up- or down-regulated between ‘US-897’ and its parents (Figure 4B). Allele-specific responses to pathogen infection were observed for three of the four quadrants for genes with a significant genotype x treatment interaction (Figure 4D; paired t-test, p < 0.01) and two of four quadrants for genes with significant treatment effect (Supplementary Figure S1D; paired t-test, p < 0.01). Overall, the signature of cis-regulation of pathogen-responsive genes is reduced compared to genes with significant genotype effects (Figure 3B), but cis-regulation still contributes to the transcriptional response to CLas infection.
Cis-regulation of genes involved in plant immunity
Genes with the largest magnitude of allele-specific expression in ‘US-897’ were additively inherited, possibly due to large-effect cis-regulatory variants. Notably, genes with additive inheritance favoring expression of the ‘Cleopatra’ allele were significantly enriched for genes involved in systemic acquired resistance (SAR) (two out of seven annotated genes; p < 0.01; Supplementary Table S6). These two genes are homologs of NONEXPRESSOR OF PATHOGENESIS-RELATED3 (AtNPR3) and NONEXPRESSOR OF PATHOGENESIS-RELATED4 (AtNPR4) in Arabidopsis thaliana and are negative regulators of immunity and critical components of plant defense against pathogens (Ding et al., 2018; Zhang et al., 2006). NPR3 has also been shown to negatively regulate immunity in sweet orange, where silencing of NPR3 expression increased basal callose levels and reduced pathogen-induced callose deposition and ROS accumulation (Sarkar et al., 2024). The expression of both NPR3 and NPR4 in ‘US-897’ was intermediate between ‘Cleopatra’ and ‘Flying Dragon’ (Figure 5A). We reasoned that reduced expression of NPR3 and NPR4 in ‘US-897’ could enhance basal immunity in the hybrid compared to its susceptible parent ‘Cleopatra’. To understand how novel combinations of cis-regulatory alleles affect expression of these genes, we evaluated allele-specific expression in ‘US-897’. Surprisingly, intermediate expression of NPR3 in ‘US-897’ results from high expression of the ‘Flying Dragon’ allele and low expression of the maternal ‘Cleopatra’ allele (Figure 5B, FDR corrected p < 0.01). This contrasts with the parental expression patterns, where ‘Cleopatra’ expression is high and ‘Flying Dragon’ expression is low (Figure 5A). The patterns of allelic expression differ at NPR4 with low expression of the ‘Flying Dragon’ allele and high expression of the ‘Cleopatra’ allele resulting in intermediate expression of the gene in ‘US-897’ compared to its two parents (Figure 5C). Patterns of allelic expression of NPR3 and NPR4 in ‘US-897’ are consistent for uninfected and CLas infected treatments. NPR3 and NPR4 negatively regulate NPR1, but because this regulation occurs post-transcriptionally (Fu et al., 2012), no relationship between levels of gene expression between NPR1 and NPR3/4 are expected. We did find that ‘Flying Dragon’ and ‘US-897’ have significantly lower NPR1 expression than ‘Cleopatra’ (Supplementary Figure S3). Overall, expression of NPR3 and NPR4 is reduced in ‘US-897’ compared to its susceptible parent, ‘Cleopatra’, and both genes have signatures of cis-regulation.

Figure 5. Cis-regulation of NPR3 and NPR4 expression in ‘US-897’. (A) The expression of homologs of NPR3 (left) and NPR4 (right) in ‘Flying Dragon’ (n=4), ‘US-897’ (n=8), and ‘Cleopatra’ (n=8). These genes are significantly differentially expressed in ‘US-897’ compared to both of its parents (FDR adjusted p-value < 0.05). (B) Allele-specific expression of NPR3 in ‘US-897’ in mock and infected treatments. (C) Allele-specific expression of NPR4 in ‘US-897’ in mock and infected treatments.
Discussion
Parent-offspring trios are commonly used to infer patterns of inheritance and, in the case of gene expression, are essential for dissecting cis- and trans- gene regulatory architecture (Hill et al., 2021). We selected a parent-offspring trio with relevance to a current major disease impacting worldwide production of citrus (Albrecht and Bowman, 2011; Court and Rahmani, 2017). The interspecific hybrid, ‘US-897’ is a commercial cultivar with tolerance to HLB disease (Albrecht and Bowman, 2011). Analysis of differential gene expression between parents and hybrid and allele-specific expression in the hybrid was made possible by the production of a phased, diploid assembly of ‘US-897’. The combination of PacBio HiFi sequencing with kmer-based trio-binning resolved both parental haplotypes of the ‘US-897’ genome, resulting in the first published reference genome containing both Citrus and Poncirus ancestry.
Genes with additive and dominant expression are cis-regulated
Expression level dominance was the most prevalent pattern of gene expression in ‘US-897’ compared to its parents. More than 80% of genes with significant differential expression between genotypes, in both untreated (Figure 2C) and treated samples (Supplementary Figure S2), were dominantly expressed in the hybrid. This is consistent with other parent-offspring trios, including in Arabidopsis thaliana F1 hybrids, where expression dominance was also detected at the majority of genes differentially expressed between hybrid and parents (Yuan et al., 2023). Interestingly, in untreated samples, dominant gene expression was slightly biased towards expression of the maternal allele (‘Cleopatra’) (Figure 2C), but this trend reversed in CLas treated samples, instead favoring dominant expression of the paternal allele (‘Flying Dragon’) (Supplementary Figure S2).
Allele-specific expression of parental alleles in hybrid genomes is an indicator of cis-regulatory control of gene expression. Analysis of allele-specific expression in ‘US-897’ revealed that a majority of genes differentially regulated between ‘Cleopatra’ and ‘Flying Dragon’ are regulated in cis. In addition to the large number of genes expressed in ‘US-897’ with ASE (59%), there was a strong correlation between ASE in ‘US-897’ and differential expression between the two parents (Figure 3A). Together this suggests that divergence in gene expression between the Citrus reticulata (‘Cleopatra’) and Poncirus trifoliata (‘Flying Dragon’) primarily results from cis-regulatory variation. Genes with the greatest magnitude of ASE in ‘US-897’ were those genes with additive inheritance and dominant inheritance with higher expression of the ‘Flying Dragon’ allele (Figure 3B). Strong allele-specific expression of these genes could be caused by large-effect mutations in cis-regulatory regions, and mining the regulatory regions of these genes could reveal regulatory motifs that have diverged between Citrus and Poncirus.
We focused on additively inherited genes with strong allele specific expression and discovered significant enrichment for genes involved in systemic acquired resistance (SAR), including homologs of Arabidopsis thaliana NONEXPRESSOR OF PATHOGENESIS-RELATED3 (AtNPR3) and NONEXPRESSOR OF PATHOGENESIS-RELATED3 (AtNPR4). These genes have established roles as negative regulators of plant immune responses (Ding et al., 2018; Zhang et al., 2006), including recently in citrus where silencing of NPR3 increased tolerance to HLB in sweet orange, suggesting that modulating NPR3 activity could tune the immune response to mitigate disease (Sarkar et al., 2024). Reduced expression of NPR3 and NPR4 in ‘US-897’ compared to susceptible ‘Cleopatra’ suggests that the hybrid may contain a favorable immune status compared to its susceptible parent. Further dissecting allele-specific expression (ASE) of NPR3 and NPR4 in ‘US-897’ revealed complex underlying regulatory architectures. For NPR4, we found evidence of a classic example of cis-regulation, with low-expression of the ‘Flying Dragon’ allele and higher expression of the ‘Cleopatra’ allele resulting in intermediate expression of NPR4 in the hybrid. In contrast, the regulatory control of NPR3 was more complex with intermediate expression of NPR3 in ‘US-897’ resulting from unexpectedly high expression of the ‘Flying Dragon’ and low expression of the ‘Cleopatra’ allele. This suggests an interplay of cis- and trans-acting factors, but the expression of each allele would need to be evaluated in both parental backgrounds to further dissect the regulation of NPR3. The precise implications and physiological impact of cis-regulation of NPR3 and NPR4 and their effect on HLB tolerance will require further functional evaluation.
Transgressive gene expression was limited in ‘US-897’
Genes with transgressive expression should not exhibit ASE, because trans-regulators will affect the expression of both parental alleles in the hybrid. Indeed, genes with transgressive expression in ‘US-897’ were not expressed in an allele-specific manner (Figure 3B). Overall, very few genes were transgressively expressed in both untreated and CLas-treated samples (Figure 2C; Supplementary Figure S2). Transgressive expression could be a consequence of hemizygosity, with the hybrid inheriting two or zero copies of a gene and resulting in transgressive up or downregulation in ‘US-897’, respectively. In fact, 10 - 15% of genes are hemizygous in the clonally propagated grapevine (Vitis vinifera) and cassava (Manihot esculenta) (Zhou et al., 2019; Peng et al., 2025). Similarly, in ‘US-897’, more than 6,000 genes were present in only one of the two parental haplotypes, making this a promising avenue for future research. Genes with transgressive expression have long been of interest for their potential contribution to hybrid vigor, but overall we find that transgressive expression is limited in this interspecific hybrid.
Genes responsive to CLas infection have a reduced signature of cis-regulation
Cis-regulatory sequences are also important for regulating gene expression in response to environmental signals (Reynoso et al., 2019; Zeng et al., 2019; Zhou et al., 2021). For example, in two species of Arabidopsis more than 6,000 drought responsive cis-regulatory variants were identified (He et al., 2021). To explore the contribution of cis-regulation to the response to CLas infection, genes with a significant treatment effect (n=629) or genotype x treatment interaction (n=1,071) were identified. Compared to previous studies of gene expression in ‘US-897’ and its maternal parent 7–8 months after infection (Albrecht and Bowman, 2012b), relatively few CLas responsive genes were identified at just 18 days after bud initiation.
Cis-regulation of genes responsive to CLas infection was limited compared to cis-regulation of gene expression divergence between ‘Cleopatra’ and ‘Flying Dragon’. This is supported by the reduced correlation of ΔASE in ‘US-897’ (the allelic response to CLas infection) to parental responses to CLas infection (Figure 4C; Supplementary Figure S1C). Also, evaluation of allelic response to CLas revealed that there were significant differences between the response of the maternal or paternal allele in ‘US-897’, although the differences were subtle (Figure 4D). This suggests that cis-regulatory variants are important for the transcriptional response to pathogen infection, but their contribution may be reduced for pathogen-responsive genes.
Relevance of genetic architecture for breeding new citrus cultivars
Many citrus rootstocks in current use for citriculture throughout the world are hybrids of Citrus species with Poncirus trifoliata (Bowman and Joubert, 2020). Poncirus trifoliata has been highly valued as a parent to create new rootstocks because of several important traits that are transmitted to progeny, including nucellar polyembryony (Kepiro and Roose, 2009), tolerance of citrus nematodes (Ling et al., 2000; Imran Hamid et al., 2025), tolerance of Phytophthora root rot (Graham, 1995), resistance to citrus tristeza virus (Fang et al., 1998), and cold tolerance (Yelenosky, 2011). More recently, P. trifoliata has also been recognized as transmitting CLas-resistance and tolerance to progeny that are useful as rootstocks (Albrecht and Bowman, 2011, 2012b), and breeding work has expanded to include P. trifoliata ancestry in advanced-generation rootstock hybrids with tolerance to CLas (Bowman et al., 2021; Bowman, 2023).
Although less widely recognized, P. trifoliata has also been incorporated in citrus scion breeding for decades, to introduce resistance to tristeza virus and better cold hardiness into new scions (Young et al., 1982; Tignor et al., 1998). More recently, the aim of such efforts has shifted to include breeding new scion hybrids that contain tolerance to CLas from P. trifoliata ancestry (Huang et al., 2023a). In citrus scion breeding, first generation hybrids with trifoliate orange typically have significant undesirable flavor components in fruit juice (Deterre et al., 2021). However, backcrossing of these P. trifoliata orange hybrids to citrus for one or more generations, combined with progeny selection, can eliminate these off-flavors, and recover new hybrids that combine the desired CLas-tolerance traits with good-flavored fruit (Fang et al., 1998; Jeffries et al., 2024). In breeding of both citrus rootstocks and citrus scions, information about regulation of gene expression for critical traits is severely limited, and lack of knowledge about transmission of key traits from parents to progeny greatly hampers progress in developing new cultivars with those essential traits. A clear understanding about genetic architecture and the regulation of gene expression associated with CLas infection will dramatically improve opportunities to make meaningful gains in the breeding of both rootstocks and scions with tolerance to CLas. While our transcriptomic analyses provide insights into the regulation of gene expression in response to CLas infection, the molecular basis of HLB tolerance is still unknown and will require validation of the roles of specific regulatory variants and genes in conferring HLB tolerance.
Conclusion
Cis-regulatory variants are a common and significant contributor to expression divergence between Citrus and its wild relative Poncirus. These cis-acting variants are associated with additive and dominant gene expression inheritance. But, genes responding to CLas shortly after infection have a reduced signature of cis-regulation, suggesting a dynamic interplay between cis- and trans-regulation of response to CLas infection in the interspecific hybrid, ‘US-897’. The high-quality, haplotype-phased genome assembly of ‘US-897’, alongside the insights gained from analysis of allele-specific expression in the ‘US-897’, provide a critical resource for uncovering the genetic basis of HLB tolerance in Citrus x Poncirus hybrids.
Materials and methods
Plant cultivation and CLas inoculation
Candidatus liberibacter asiaticus (CLas) is a quarantined pathogen in California and all plant cultivation and psyllid-mediated inoculations were performed in the California Citrus Research Foundation containment facility located in Riverside, CA. Seeds of ‘Cleopatra’ mandarin (Citrus reticulata), ‘Flying Dragon’ trifoliate orange (Poncirus trifoliata), and their F1 hybrid ‘US-897’ were sterilized prior to germination. The primary seed coats were removed, and seeds were washed with 70% isopropanol for 10 minutes. Next, seeds were treated with a 10% bleach solution for 30 minutes and rinsed three times with sterile water. Seeds were placed in 15 ml culture tubes containing seed germination media (0.22% Murashige and Skoog basal salt mixture (Murashige and Skoog, 1962), 2.5% sucrose, 0.005% myo-inositol, 0.50% ethylenediaminetetraacetic acid ferric salt, and 0.70% agar) and grown at 30° C in the dark for two weeks. After germination, seeds were transplanted to cone-tainers (164 ml) with coconut coir and placed in a greenhouse with 12 hour supplemental light:dark cycles with light intensity maintained at ~ 500 µmol/m²/s. Both ‘Cleopatra’ and ‘Flying Dragon’ produce nucellar seeds. All transplanted seedlings were genotyped with 24 LGC KASP™ markers (Hoddesdon, UK) and zygotic seedlings were culled.
For CLas inoculation, Asian citrus psyllids (Diaphorina citri) were reared on either confirmed CLas-infected citron (Citrus medica) plants or CLas-free citron plants, thereby providing CLas-positive (CLas+) and mock-inoculation (CLas-) psyllid populations, respectively. Nine month old plants (6 replicates) of each genotype (‘Cleopatra’, ‘Flying Dragon’, US-897) were infested with 15–20 CLas+ psyllids per plant for a period of four weeks. Concurrently, three replicate plants of each genotype were similarly infested with 15–20 clean (CLas-) psyllids per plant to serve as mock-inoculated controls. Psyllids were contained on the plants using small mesh cages enclosing a single branch during the infestation period. In some rare cases, smaller plants were completely caged. Plants were nine months-old at the time of infestation with psyllids. Before tissue sampling for RNA-sequencing, all experimental plants were pruned to stimulate new flush. At this time, the systemic CLas infection status of each individual plant was confirmed by quantitative PCR of the nrdB gene of CLas as performed in (Zheng et al., 2016). The leaf sampled to confirm CLas infection was from the single, caged branch (or small plant). This branch was trimmed and then bud initiation was tracked, with sampling of leaf tissue for RNA-seq occurring 18 days after bud initiation. Based on qPCR results, five ‘Cleopatra’, two ‘Flying Dragon’, and five US-897 individuals from the CLas+ treatment group were confirmed as infected and used for subsequent RNA-seq analysis (Supplementary Table S2). All plants, mock- and CLas-infected, were randomized in both the plant growth chamber, where psyllid inoculations occurred, and in the greenhouse.
Sampling scheme and RNA-sequencing
Plants were transferred from the growth chamber used for psyllid-mediated CLas (and mock) inoculations to the greenhouse two months after their initial exposure to psyllids. Plants continued to grow for an additional two months in the greenhouse and then, after adjusting to greenhouse conditions, leaf samples were collected to assay for CLas. This enables sufficient establishment of the bacteria in the host plants. Additionally, CLas detection in new shoots has been shown to occur as early as 15 days post-bud initiation (Ma et al., 2022). After sampling, a single branch on each plant was trimmed to induce new flush and plants were monitored daily for bud initiation. Next, leaf tissue from the newly emerging flush was sampled 18 days post bud initiation to capture early responses to pathogen infection as performed in (Ma et al., 2022). While the exact newly emerging shoots sampled for RNA-seq were not individually quantified for bacterial titer, CLas is known to be transported to and can be detected in newly developing flush tissues within 15 days of bud-initiation (Ma et al., 2022). Therefore, RNA-seq was performed on these newly emerged shoots from systemically infected plants. All sampling was performed from 10-11AM Pacific Standard Time (PST). After sampling, leaf tissue was flash frozen in liquid nitrogen and quarantined at -80° C for 14 days. RNA extraction was then performed using the Qiagen RNeasy Plant kit (Hilden, Germany). Sequencing libraries were prepared using the Illumina Stranded mRNA prep, Ligation kit (Pleasanton, California). Libraries were sequenced on two lanes of a NovaSeq X Plus (PE150) by Novogene (Beijing, China) generating an average of 28.3 million reads per sample (Supplementary Tables S3, S4).
PacBio sequencing
Nuclei extraction from ~3 grams of young leaf tissue from a single individual of the F1 hybrid ‘US-897’ was performed using PacBio’s nuclei extraction protocol (PacBio, 2022b). High molecular weight DNA (HMW DNA) was then extracted from nuclei using PacBio’s Nanobind plant kit (PacBio, 2022a). An extended lysis step was included, increasing from 30 minutes to 2 hours. A second lysis step was performed to ensure sufficient extraction of HMW DNA. HMW DNA was then sent to Novogene (Beijing, China) for library preparation and sequencing on one SMRT cell of the PacBio Revio (Menlo Park, California). This generated 65 Gb of HiFi sequencing reads, with an average read length of 15.6 Kb.
Genome-assembly and annotation
PacBio HiFi sequencing reads were used to assemble the diploid genome of the F1 hybrid ‘US-897’. First, Illumina whole-genome sequencing reads from the parental genotypes, ‘Cleopatra’ mandarin and ‘Flying Dragon’ trifoliate orange (Wu et al., 2018; Peng et al., 2020), were used to identify parent-specific k-mers (k=19), or hap-mers, using Merqury (Rhie et al., 2020). PacBio HiFi reads were error-corrected and used to construct an assembly graph using verkko (Rautiainen et al., 2023). Next, paths in this assembly graph were haplotype-resolved using parental hap-mer profiles (verkko trio), and a consensus algorithm was applied to construct the final diploid assembly (Rautiainen et al., 2023). Before genome annotation, repetitive sequences were soft-masked using a custom library of repeat sequences with RepeatMasker v 4.1.5 (Smit et al., 2013-2015). Genome annotation was then performed with BRAKER v3.0.3, using a custom database of proteins, including those specific to the genus Citrus. Only the 9 largest scaffolds of each haplotype were annotated. BUSCO scores were calculated for the genome assembly (-m geno) and the genome annotation (-m prot) using the eudicots_db10 database and Busco v5.8.0 (Manni et al., 2021).
Differential gene expression analysis
Genome annotations for each haplotype of ‘US-897’ were used to construct a diploid transcriptome using Salmon v1.7.0 (Patro et al., 2017). Homology of transcripts between haplotypes was inferred based on protein sequence identity and collinearity using GeneTribe (Chen et al., 2020). Transcripts were classified based on whether they are reciprocal best hits (RBH) that are allelic between haplotypes, single best hits (SBH) are orthologous sequences, and singletons are genes that are unique to a haplotype. For ortholog identification, the gene IDs from the ‘Flying Dragon’ haplotype were used, unless the gene was only present in Cleopatra, then the gene IDs from the ‘Cleopatra’ haplotype were used. Next transcripts were quantified using Salmon v1.7.0 (Patro et al., 2017) with the following parameters: -validateMappings –numBootstraps 30 –recoverOrphans –writeUnmappedNames –allowDovetail –softclip –minScoreFraction 0.6. Bootstrapping (n=30) was performed to generate inferential replicates for accurate quantification of transcripts from RNA-sequencing reads. For analyses of overall gene-expression (as opposed to allele-specific expression), read counts for transcripts were summed for the set of orthologous genes. Gene expression estimates (Counts-Per-Million or CPM) were log2-transformed and weights were assigned to each observation based on the estimated mean-variance trend across genes using voom transformation (Law et al., 2014). Next, gene expression estimates for the estimated mid-parent values for each gene were calculated by averaging read counts of parental samples that were randomly selected from replicates of each treatment. This process was repeated to generate four replicates of mid-parent expression per gene. Genes with low levels of expression were removed (CPM less than 2) and the remaining set of 20,981 genes were the basis of differential gene expression analysis. Gene expression estimates were modeled using limma-voom (Law et al., 2014), including fixed effects for genotype (‘Cleopatra’, ‘Flying Dragon’, ‘US-897’, mid-parent), treatment (CLas-inoculated, mock-inoculated), and the genotype x treatment interaction. This mid-parent value-based classification, performed using limma’s linear models and statistical comparisons (Ritchie et al., 2015), serves as a standard approach for broadly categorizing gene expression inheritance patterns (Rapp et al., 2009; Huang et al., 2016). While the midparent value establishes a null hypothesis for additive inheritance for each gene, the use of limma accounts for intra-genotypic variability in its variance estimation, thereby providing robust statistical significance for deviations.
Classification of hybrid gene expression inheritance patterns
Specific contrasts were evaluated to identify differentially expressed genes (DEGs) between genotypes (‘Cleopatra’ versus ‘Flying Dragon’, ‘Cleopatra’ or ‘Flying Dragon’ versus ‘US-897’, and ‘US-897’ versus estimated mid-parent value). Gene expression inheritance patterns were determined for genes with differential expression in any of these contrasts using the R package “HybridExpress” (Almeida-Silva et al., 2023). Significant differences based on Benjamini-Hochberg FDR < 0.05 and the log2(fold-change) for each contrast were used to classify gene-expression patterns. Genes were classified as having “additive” inheritance if there was no significant difference in expression (Benjamini-Hochberg FDR > 0.05) between the expected mid-parent value and ‘US-897’. Genes were classified as “dominant” if there was a significant difference between the hybrid and one parent, but not the other. Lastly, genes with significantly higher or lower gene expression in the ‘US-897’ compared to both parents were classified as transgressive. This same classification scheme was used to assess gene expression inheritance patterns for genotype x treatment contrasts.
Identification of allele-specific expression in ‘US-897’
Gene expression was quantified as described above using Salmon, with bootstrapping (n=30) to generate inferential replicates for accurate quantification of alleles from RNA-sequencing reads. Gene expression was quantified for 22,419 genes, and genes with read counts greater than ten in at least six samples were retained (n=17,487). Genes for which inferential replicates of each allele were almost identical were filtered from downstream analyses (n=2), as this indicates there is no information about allelic expression in the sequencing reads. For each gene, the ratio of allele 1 (‘Flying Dragon’): allele 2 (‘Cleopatra’) expression was averaged across all samples, and genes with a ratio greater than 0.95 (n=371) or less than 0.05 (n=467) were classified as having monoallelic expression. Monoallelic genes were excluded from modeling of allele-specific expression resulting in 15,549 genes that could be evaluated for allele-specific expression.
Data availability statement
All raw sequencing data has been deposited at NCBI, including RNA-sequencing reads (BioSamples: SAMN50643246-SAMN50643264, SAMN506432650) and PacBio HIFI reads used to assemble the genome of 'US-897' (BioSample: SAMN50643266). The genome assembly of 'US-897' has also been deposited at NCBI under BioProject PRJNA1293460 and PRJNA1293461 for the paternal and maternal haplotypes, respectively.
Author contributions
ID: Conceptualization, Data curation, Formal Analysis, Investigation, Methodology, Validation, Visualization, Writing – original draft, Writing – review & editing. OZ: Investigation, Methodology, Writing – review & editing. EÁ: Formal Analysis, Investigation, Methodology, Writing – review & editing. CJ: Investigation, Methodology, Writing – review & editing. KB: Conceptualization, Funding acquisition, Investigation, Resources, Writing – review & editing. DS: Conceptualization, Funding acquisition, Investigation, Project administration, Resources, Supervision, Writing – original draft, Writing – review & editing.
Funding
The author(s) declare financial support was received for the research and/or publication of this article. This research was supported by USDA-NIFA Award # 2020-70029-33202 to DS; USDA-NIFA Award # 2023-70029-41305 to DS and KB; ID is a fellow in the Plants-3D NSF National Research Traineeship Program (DBI-1922642). Computations were performed using the computer clusters and data storage resources of the UC Riverside HPCC, which were funded by grants from NSF (MRI-2215705, MRI-1429826) and NIH (1S10OD016290-01A1).
Acknowledgments
We thank Dr. Le’Kneitah Smith and research staff at the California Citrus Research Foundation Contained Research Facility for their guidance and support of research on Huanglongbing disease in citrus.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that no Generative AI was used in the creation of this manuscript.
Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2025.1627531/full#supplementary-material
Supplementary Figure 1 | Classification of patterns of gene expression inheritance in ‘Cleopatra’, ‘US-897’, and ‘Flying Dragon’ for genes with significant treatment effects (n=629). (A) Identification of genes with significant treatment effects. The dotted line represents the significance threshold (FDR corrected p < 0.05). (B) Categorization of gene expression patterns relative to expression in ‘US-897’ versus parental genotypes for genes with significant genotype x treatment interaction (n=629). The dot-plot compares the log2 fold-change in gene expression between infected samples of ‘US-897’ and its two parents, with points colored by quadrant. (C) The relationship between Δ ASE in US-897 (the allelic response to CLas infection) (log2(FD/CLEO)Allelic) and expression divergence between parental genotypes ‘Flying Dragon’ and ‘Cleopatra’ (log2(FD/CLEO)Parental). Only genes with significant treatment effects that could be tested for ASE are included (n=535). (D) The response of each allele in ‘US-897’ to pathogen infection for genes in each of the quadrants in (B). Only genes with significant treatment effects that could be tested for ASE are included (n=535).
Supplementary Figure 2 | Classification of patterns of gene expression inheritance in ‘Cleopatra’, ‘US-897’, and ‘Flying Dragon’ in untreated (left, n=5,949) and treated (right, n=5,856) samples. (A) Diagram of the number of differentially expressed genes (DEGs) with a significant genotype effect in comparisons of: ‘Cleopatra’ versus ‘Flying Dragon’, ‘Cleopatra’ or ‘Flying Dragon’ versus ‘US-897’, and ‘US-897’ versus estimated mid-parent value. Solid black lines indicate comparison, with the total number of DEGs per comparison noted in the center of each line and the number of upregulated DEGs per comparison noted adjacent to each genotype. Percentages are based on the total number of tested genes (n = 20,981). (B) The 12 categories of gene expression inheritance were grouped into five major classifications: transgressive up in ‘US-897’, transgressive down in US-897, additive, dominant for ‘Flying Dragon’ allele, and dominant for ‘Cleopatra’ allele.
Supplementary Figure 3 | NPR1 expression in ‘US-897’ and its parents. The expression of homologs of NPR1 in ‘Flying Dragon’ (n=4), ‘US-897’ (n=8), and ‘Cleopatra’ (n=8) under mock and CLas. NPR1 is significantly differentially expressed in ‘US-897’ compared to ‘Cleopatra’ (FDR adjusted p-value < 0.05) but not ‘Flying Dragon’.
References
Albrecht, U. and Bowman, K. D. (2011). Tolerance of the Trifoliate Citrus Hybrid US-897 (Citrus reticulata Blanco × Poncirus trifoliata L. Raf.) to Huanglongbing. HortScience 46, 16–22. doi: 10.21273/HORTSCI.46.1.16
Albrecht, U. and Bowman, K. D. (2012a). Tolerance of trifoliate citrus rootstock hybrids to Candidatus Liberibacter asiaticus. Sci. Hortic. 147, 71–80. doi: 10.1016/j.scienta.2012.08.036
Albrecht, U. and Bowman, K. D. (2012b). Transcriptional response of susceptible and tolerant citrus to infection with Candidatus Liberibacter asiaticus. Plant Sci. 185-186, 118–130. doi: 10.1016/j.plantsci.2011.09.008
Albrecht, U., Tripathi, I., and Bowman, K. D. (2020). Rootstock influences the metabolic response to Candidatus Liberibacter asiaticus in grafted sweet orange trees. Trees 34, 405–431. doi: 10.1007/s00468-019-01925-3
Almeida-Silva, F., Prost-Boxoen, L., and Van de Peer, Y. (2023). HybridExpress: Comparative analysis of RNA-seq data for hybrids and their progenitors. doi: 10.1111/nph.19862
Bové, J. M. (2006). Huanglongbing: A destructive, newly-emerging, century-old disease of citrus. J. Plant Pathol. 88, 7–37.
Bowman, K. D. (2023). Four new SuperSour rootstocks for improved production of sweet orange in a huanglongbing environment. HortScience 58, 1622–1625. doi: 10.21273/HORTSCI17466-23
Bowman, K. D. and Albrecht, U. (2020). Rootstock influences on health and growth following candidatus liberibacter asiaticus infection in young sweet orange trees. Agronomy 10, 1907. doi: 10.3390/agronomy10121907
Bowman, K. D. and Joubert, J. (2020). “Citrus rootstocks,” in The Genus Citrus (Cambridge, Mass: Elsevier), 105–127.
Bowman, K. D. and McCollum, G. (2015). Five new citrus rootstocks with improved tolerance to huanglongbing. HortScience 50, 1731–1734. doi: 10.21273/HORTSCI.50.11.1731
Bowman, K. D., McCollum, G., and Albrecht, U. (2021). SuperSour: A new strategy for breeding superior citrus rootstocks. Front. Plant Sci. 12, 741009. doi: 10.3389/fpls.2021.741009
Chen, Z. J. (2013). Genomic and epigenetic insights into the molecular bases of heterosis. Nat. Rev. Genet. 14, 471–482. doi: 10.1038/nrg3503
Chen, Y., Song, W., Xie, X., Wang, Z., Guan, P., Peng, H., et al. (2020). A collinearity-incorporating homology inference strategy for connecting emerging assemblies in the triticeae tribe as a pilot practice in the plant pangenomic era. Mol. Plant 13, 1694–1708. doi: 10.1016/j.molp.2020.09.019
Coolon, J. D., McManus, C. J., Stevenson, K. R., Graveley, B. R., and Wittkopp, P. J. (2014). Tempo and mode of regulatory evolution in Drosophila. Genome Res. 24, 797–808. doi: 10.1101/gr.163014.113
Court, C. D. and Rahmani, E. S. M. (2017). Economic Contributions of the Florida Citrus Industry in 2015-2016. Available online at: https://fred.ifas.ufl.edu/media/fredifasufledu/photos/economic-impact/Economic-Impacts-of-the-Florida-Citrus-Industry-2015-16.pdf (Accessed January 31, 2025).
Deterre, S. C., McCollum, G., Leclair, C., Manthey, J. A., Bai, J., Baldwin, E. A., et al. (2021). Effect of Poncirus trifoliata on the chemical composition of fruits in pedigrees of Citrus scion hybrids. Sci. Hortic. (Amsterdam) 277, 109816. doi: 10.1016/j.scienta.2020.109816
Ding, Y., Sun, T., Ao, K., Peng, Y., Zhang, Y., Li, X., et al. (2018). Opposite roles of salicylic acid receptors NPR1 and NPR3/NPR4 in transcriptional regulation of plant immunity. Cell 173 (6), 1454-1467. doi: 10.1016/j.cell.2018.03.044
Emerson, J. J., Hsieh, L.-C., Sung, H.-M., Wang, T.-Y., Huang, C.-J., Lu, H. H.-S., et al. (2010). Natural selection on cis and trans regulation in yeasts. Genome Res. 20, 826–836. doi: 10.1101/gr.101576.109
Fang, D. Q., Federici, C. T., and Roose, M. L. (1998). A high-resolution linkage map of the citrus tristeza virus resistance gene region in Poncirus trifoliata (L.) Raf. Genetics 150, 883–890. doi: 10.1093/genetics/150.2.883
Fu, Z. Q., Yan, S., Saleh, A., Wang, W., Ruble, J., Oka, N., et al. (2012). NPR3 and NPR4 are receptors for the immune signal salicylic acid in plants. Nature 486, 228–232. doi: 10.1038/nature11162
Garcia, A. A. F., Wang, S., Melchinger, A. E., and Zeng, Z.-B. (2008). Quantitative trait loci mapping and the genetic basis of heterosis in maize and rice. Genetics 180, 1707–1724. doi: 10.1534/genetics.107.082867
Goncalves, A., Leigh-Brown, S., Thybert, D., Stefflova, K., Turro, E., Flicek, P., et al. (2012). Extensive compensatory cis-trans regulation in the evolution of mouse gene expression. Genome Res. 22, 2376–2384. doi: 10.1101/gr.142281.112
Graham, J. (1995). Root regeneration and tolerance of citrus rootstocks to root rot caused byPhytophthora nicotianae. Phytopathology 85, 111.
Guo, M., Rupe, M. A., Wei, J., Winkler, C., Goncalves-Butruille, M., Weers, B. P., et al. (2014). Maize ARGOS1 (ZAR1) transgenic alleles increase hybrid maize yield. J. Exp. Bot. 65, 249–260. doi: 10.1093/jxb/ert370
Hall, D. G. (2008). BIOLOGY, HISTORY AND WORLD STATUS OF Diaphorina citri. Available online at: https://www.ars.usda.gov/ARSUserFiles/35403/HallDiaphorinabiologyhistoryworldstatus2008.pdf (Accessed October 6, 2022).
Hashimoto, S., Wake, T., Nakamura, H., Minamiyama, M., Araki-Nakamura, S., Ohmae-Shinohara, K., et al. (2021). The dominance model for heterosis explains culm length genetics in a hybrid sorghum variety. Sci. Rep. 11, 4532. doi: 10.1038/s41598-021-84020-3
He, F., Steige, K. A., Kovacova, V., Göbel, U., Bouzid, M., Keightley, P. D., et al. (2021). Cis-regulatory evolution spotlights species differences in the adaptive potential of gene expression plasticity. Nat. Commun. 12, 3376. doi: 10.1038/s41467-021-23558-2
Hill, M. S., Vande Zande, P., and Wittkopp, P. J. (2021). Molecular and evolutionary processes generating variation in gene expression. Nat. Rev. Genet. 22, 203–215. doi: 10.1038/s41576-020-00304-w
Hua, J., Xing, Y., Wu, W., Xu, C., Sun, X., Yu, S., et al. (2003). Single-locus heterotic effects and dominance by dominance interactions can adequately explain the genetic basis of heterosis in an elite rice hybrid. Proc. Natl. Acad. Sci. U. S. A. 100, 2574–2579. doi: 10.1073/pnas.0437907100
Huang, Y., He, J., Xu, Y., Zheng, W., Wang, S., Chen, P., et al. (2023b). Pangenome analysis provides insight into the evolution of the orange subfamily and a key gene for citric acid accumulation in citrus fruits. Nat. Genet. 55, 1964–1975. doi: 10.1038/s41588-023-01516-6
Huang, M., Roose, M. L., Yu, Q., Stover, E., Hall, D. G., Deng, Z., et al. (2023a). Mapping of QTLs and candidate genes associated with multiple phenotypic traits for Huanglongbing tolerance in citrus. Hortic. Plant J. 9, 705–719. doi: 10.1016/j.hpj.2022.10.008
Huang, X., Yang, S., Gong, J., Zhao, Q., Feng, Q., Zhan, Q., et al. (2016). Genomic architecture of heterosis for yield traits in rice. Nature 537, 629–633. doi: 10.1038/nature19760
Imran Hamid, M., Bowman, K. D., Seymour, D. K., and Rolshausen, P. E. (2025). Evaluation of huanglongbing-tolerant citrus hybrid rootstocks for resistance to tylenchulus semipenetrans. Plant Health Progress. doi: 10.1094/PHP-08-24-0074-RS
Jeffries, K. A., Fan, Z., Sun, X., Olmedo, G. M., Zhao, W., Mattia, M., et al. (2024). New insights in the flavor and chemistry of Huanglongbing tolerant citrus hybrids with/without Poncirus trifoliata in their pedigree. Front. Hortic. 3. doi: 10.3389/fhort.2024.1425366
Kepiro, J. L. and Roose, M. L. (2009). AFLP markers closely linked to a major gene essential for nucellar embryony (apomixis) in Citrus maxima × Poncirus trifoliata. Tree Genet. Genomes 6, 1–11. doi: 10.1007/s11295-009-0223-z
Krieger, U., Lippman, Z. B., and Zamir, D. (2010). The flowering gene SINGLE FLOWER TRUSS drives heterosis for yield in tomato. Nat. Genet. 42, 459–463. doi: 10.1038/ng.550
Law, C. W., Chen, Y., Shi, W., and Smyth, G. K. (2014). voom: Precision weights unlock linear model analysis tools for RNA-seq read counts. Genome Biol. 15, R29. doi: 10.1186/gb-2014-15-2-r29
Li, L., Lu, K., Chen, Z., Mu, T., Hu, Z., and Li, X. (2008). Dominance, overdominance and epistasis condition the heterosis in two heterotic rice hybrids. Genetics 180, 1725–1742. doi: 10.1534/genetics.108.091942
Ling, P., Duncan, L. W., Deng, Z., Dunn, D., Hu, X., Huang, S., et al. (2000). Inheritance of citrus nematode resistance and its linkage with molecular markers. Theor. Appl. Genet. 100, 1010–1017. doi: 10.1007/s001220051382
Ma, W., Pang, Z., Huang, X., Xu, J., Pandey, S. S., Li, J., et al. (2022). Citrus Huanglongbing is a pathogen-triggered immune disease that can be mitigated with antioxidants and gibberellin. Nat. Commun. 13, 529. doi: 10.1038/s41467-022-28189-9
Manni, M., Berkeley, M. R., Seppey, M., Simão, F. A., and Zdobnov, E. M. (2021). BUSCO update: novel and streamlined workflows along with broader and deeper phylogenetic coverage for scoring of eukaryotic, prokaryotic, and viral genomes. Mol. Biol. Evol. 38, 4647–4654. doi: 10.1093/molbev/msab199
Metzger, B. P. H., Wittkopp, P. J., and Coolon, J. D. (2017). Evolutionary dynamics of regulatory changes underlying gene expression divergence among saccharomyces species. Genome Biol. Evol. 9, 843–854. doi: 10.1093/gbe/evx035
Murashige, T. and Skoog, F. (1962). A revised medium for rapid growth and bio assays with tobacco tissue cultures. Physiol. Plant 15, 473–497. doi: 10.1111/j.1399-3054.1962.tb08052.x
PacBio (2022a). Extracting HMW DNA from plant nuclei using Nanobind® kits. Available online at: https://www.pacb.com/wp-content/uploads/Procedure-checklist-Extracting-HMW-DNA-from-plant-nuclei-using-Nanobind-kits.pdf (Accessed January 31, 2024).
PacBio (2022b). Isolating nuclei from plant tissue using LN2 disruption. Available online at: https://www.pacb.com/wp-content/uploads/Procedure-checklist-Isolating-nuclei-from-plant-tissue-using-LN2-disruption.pdf (Accessed January 31, 2024).
Paril, J., Reif, J., Fournier-Level, A., and Pourkheirandish, M. (2024). Heterosis in crop improvement. Plant J. 117, 23–32. doi: 10.1111/tpj.v117.1
Patro, R., Duggal, G., Love, M. I., Irizarry, R. A., and Kingsford, C. (2017). Salmon provides fast and bias-aware quantification of transcript expression. Nat. Methods 14, 417–419. doi: 10.1038/nmeth.4197
Peng, Z., Bredeson, J. V., Wu, G. A., Shu, S., Rawat, N., Du, D., et al. (2020). A chromosome-scale reference genome of trifoliate orange (Poncirus trifoliata) provides insights into disease resistance, cold tolerance and genome evolution in Citrus. Plant J. 104, 1215–1232. doi: 10.1111/tpj.v104.5
Peng, Y., Wang, Y., Liu, Y., Fang, X., Cheng, L., Long, Q., et al. (2025). The genomic and epigenomic landscapes of hemizygous genes across crops with contrasting reproductive systems. Proc. Natl. Acad. Sci. U. S. A. 122, e2422487122. doi: 10.1073/pnas.2422487122
Pereira, V. and Kuzmin, E. (2024). Trans-regulatory variant network contributes to missing heritability. Cell Genom. 4, 100470. doi: 10.1016/j.xgen.2023.100470
Rapp, R. A., Udall, J. A., and Wendel, J. F. (2009). Genomic expression dominance in allopolyploids. BMC Biol. 7, 18. doi: 10.1186/1741-7007-7-18
Rautiainen, M., Nurk, S., Walenz, B. P., Logsdon, G. A., Porubsky, D., Rhie, A., et al. (2023). Telomere-to-telomere assembly of diploid chromosomes with Verkko. Nat. Biotechnol. 41, 1474–1482. doi: 10.1038/s41587-023-01662-6
Reynoso, M. A., Kajala, K., Bajic, M., West, D. A., Pauluzzi, G., Yao, A. I., et al. (2019). Evolutionary flexibility in flooding response circuitry in angiosperms. Science 365, 1291–1295. doi: 10.1126/science.aax8862
Rhie, A., Walenz, B. P., Koren, S., and Phillippy, A. M. (2020). Merqury: reference-free quality, completeness, and phasing assessment for genome assemblies. Genome Biol. 21, 245. doi: 10.1186/s13059-020-02134-9
Rockman, M. V. and Kruglyak, L. (2006). Genetics of global gene expression. Nat. Rev. Genet. 7, 862–872. doi: 10.1038/nrg1964
Sarkar, P., El-Mohtar, C., Turner, D., Welker, S., Robertson, C. J., Orbovic, V., et al (2024). NONEXPRESSOR OF PATHOGENESIS-RELATED GENES control Huanglongbing tolerance by regulating immune balance in citrus plants. bioRxiv. doi: 10.1101/2024.03.18.585579v1
Schwartz, D. and Laughner, W. J. (1969). A molecular basis for heterosis. Science 166, 626–627. doi: 10.1126/science.166.3905.626
Semel, Y., Nissenbaum, J., Menda, N., Zinder, M., Krieger, U., Issman, N., et al. (2006). Overdominant quantitative trait loci for yield and fitness in tomato. Proc. Natl. Acad. Sci. U. S. A. 103, 12981–12986. doi: 10.1073/pnas.0604635103
Seymour, D. K., Chae, E., Grimm, D. G., Martín Pizarro, C., Habring-Müller, A., Vasseur, F., et al. (2016). Genetic architecture of nonadditive inheritance in Arabidopsis thaliana hybrids. Proc. Natl. Acad. Sci. U. S. A. 113, E7317–E7326. doi: 10.1073/pnas.1615268113
Signor, S. A. and Nuzhdin, S. V. (2018). The Evolution of Gene Expression in cis and trans. Trends Genet. 34, 532–544. doi: 10.1016/j.tig.2018.03.007
Smit, A. F. A., Hubley, R., and Green, P. (2013-2015). RepeatMasker. Open-4.0. Available online at: http://www.repeatmasker.org (Accessed January 31, 2024).
Snoek, B. L., Sterken, M. G., Bevers, R. P. J., Volkers, R. J. M., Van’t Hof, A., Brenchley, R., et al. (2017). Contribution of trans regulatory eQTL to cryptic genetic variation in C. elegans. BMC Genomics 18, 500. doi: 10.1186/s12864-017-3899-8
Tian, S., Xu, X., Zhu, X., Wang, F., Song, X., and Zhang, T. (2019). Overdominance is the major genetic basis of lint yield heterosis in interspecific hybrids between G. hirsutum and G. barbadense. Heredity (Edinb.) 123, 384–394. doi: 10.1038/s41437-019-0211-5
Tignor, M. E., Davies, F. S., and Sherman, W. B. (1998). Freezing tolerance and growth characteristics of USDA intergeneric citrus hybrids US 119 and selection 17-11. HortScience 33, 744–748. doi: 10.21273/HORTSCI.33.4.744
Videvall, E., Sletvold, N., Hagenblad, J., Ågren, J., and Hansson, B. (2016). Strong maternal effects on gene expression in Arabidopsis lyrata hybrids. Mol. Biol. Evol. 33, 984–994. doi: 10.1093/molbev/msv342
Vrebalov, J., Ruezinsky, D., Padmanabhan, V., White, R., Medrano, D., Drake, R., et al. (2002). A MADS-box gene necessary for fruit ripening at the tomato ripening-inhibitor (rin) locus. Science 296, 343–346. doi: 10.1126/science.1068181
Wittkopp, P. J., Haerum, B. K., and Clark, A. G. (2004). Evolutionary changes in cis and trans gene regulation. Nature 430, 85–88. doi: 10.1038/nature02698
Wu, E. Y., Singh, N. P., Choi, K., Zakeri, M., Vincent, M., Churchill, G. A., et al. (2023). SEESAW: detecting isoform-level allelic imbalance accounting for inferential uncertainty. Genome Biol. 24, 165. doi: 10.1186/s13059-023-03003-x
Wu, G. A., Terol, J., Ibanez, V., López-García, A., Pérez-Román, E., Borredá, C., et al. (2018). Genomics of the origin and evolution of Citrus. Nature 554, 311–316. doi: 10.1038/nature25447
Xiao, Y., Jiang, S., Cheng, Q., Wang, X., Yan, J., Zhang, R., et al. (2021). The genetic mechanism of heterosis utilization in maize improvement. Genome Biol. 22, 148. doi: 10.1186/s13059-021-02370-7
Yelenosky, G. (2011). “Cold Hardiness in Citrus,” in Horticultural Reviews (John Wiley & Sons, Inc, Hoboken, NJ, USA), 201–238.
Young, R., Barrett, H. C., Hearn, C. J., and Hutchison, D. J. (1982). New sources of cold hardiness for citrus breeding. HortScience 17, 866. doi: 10.21273/HORTSCI.17.6.866
Yuan, W., Beitel, F., Srikant, T., Bezrukov, I., Schäfer, S., Kraft, R., et al. (2023). Pervasive under-dominance in gene expression underlying emergent growth trajectories in Arabidopsis thaliana hybrids. Genome Biol. 24, 200. doi: 10.1186/s13059-023-03043-3
Zeng, Z., Zhang, W., Marand, A. P., Zhu, B., Buell, C. R., and Jiang, J. (2019). Cold stress induces enhanced chromatin accessibility and bivalent histone modifications H3K4me3 and H3K27me3 of active genes in potato. Genome Biol. 20, 123. doi: 10.1186/s13059-019-1731-2
Zhang, X. and Emerson, J. J. (2019). Inferring the genetic architecture of expression variation from replicated high throughput allele-specific expression experiments. Cold Spring Harbor Lab., 699074. doi: 10.1101/699074
Zhang, Y., Cheng, Y. T., Qu Zhao, N. Q., Bi, D., and Li, X. (2006). Negative regulation of defense responses in Arabidopsis by two NPR1 paralogs. Plant J. doi: 10.1111/j.1365-313X.2006.02903.x
Zheng, Z., Xu, M., Bao, M., Wu, F., Chen, J., and Deng, X. (2016). Unusual Five Copies and Dual Forms of nrdB in “Candidatus Liberibacter asiaticus”: Biological Implications and PCR Detection Application. Sci. Rep. 6, 39020. doi: 10.1038/srep39020
Zhou, P., Enders, T. A., Myers, Z. A., Magnusson, E., Crisp, P. A., Noshay, J., et al. (2021). Prediction of conserved and variable heat and cold stress response in maize using cis-regulatory information. Plant Cell 34, 514–534. doi: 10.1093/plcell/koab267
Keywords: Citrus, Poncirus, Huanglongbing, allele-specific expression (ASE), cis-regulation
Citation: Diaz IA, Zayed O, Ávila De Dios E, Jiang C, Bowman KD and Seymour DK (2025) The genetic architecture of gene expression regulation in a Citrus x Poncirus hybrid tolerant to Huanglongbing. Front. Plant Sci. 16:1627531. doi: 10.3389/fpls.2025.1627531
Received: 12 May 2025; Accepted: 31 July 2025;
Published: 04 September 2025.
Edited by:
Eduardo Augusto Girardi, Brazilian Agricultural Research Corporation (EMBRAPA), BrazilReviewed by:
Nian Wang, University of Florida, United StatesAbelmon Gesteira, Brazilian Agricultural Research Corporation (EMBRAPA), Brazil
Qibin Yu, University of Florida, United States
Copyright © 2025 Diaz, Zayed, Ávila De Dios, Jiang, Bowman and Seymour. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Danelle K. Seymour, ZHNleW1vdXJAdWNyLmVkdQ==