Long-Read Transcriptomes Reveal Purifying Selection and Species-Linked Expression Shifts in Five Rosaceae Species

A team of researchers has generated long-read reference transcriptomes for five species of the Rosaceae tribe Spiraeeae, a poorly characterized clade that includes the genera Aruncus, Spiraea, and Petrophytum. The study, published in Plant Direct, combined full-length Iso-Seq sequencing with 3′-tag RNA sequencing (TagSeq) to build a comparative resource for examining sequence divergence, species-associated expression variation, and candidate genes involved in leaf development. The work addresses a notable gap: while economically significant Rosaceae members such as apples, strawberries, and pears have been extensively studied genetically, many wild or less commercially prominent species in the family lack in-depth genetic data.

Tribe Spiraeeae, within the subfamily Amygdaloideae, is a well-supported clade comprising nine genera and 65 to 100 species that collectively display remarkable morphological and ecological diversity. Previous work by Potter and colleagues in 2007 divided the tribe into two subclades, one comprising Aruncus, Luetkea, Holodiscus, and Xerospiraea, and the other including the remaining genera. The five species sampled in the new study span a broad ecological gradient within the tribe. Petrophytum caespitosum is associated with dry, exposed rocky substrates, whereas Aruncus dioicus typically occupies mesic habitats such as forest margins and openings. Within Spiraea, S. lucida is found in relatively drier to mesic habitats, whereas S. splendens is associated with wetter habitats, including streamsides and marshy alpine meadows. The leaf morphologies of the sampled species range from large and highly compound with toothed leaflet blades in Aruncus, to moderately small and simple leaves with toothed margins in Spiraea, to very small and entire leaves in Petrophytum.

Plant materials included three to six biological replicates of each species, harvested either from plants growing in the field for P. caespitosum subsp. caespitosum and S. splendens or from potted plants for Aruncus aethusifolius, A. dioicus, and S. lucida. Field collections were made in Salmon–Challis National Forest in Idaho and Lassen National Forest in California. The three to five youngest visible leaves at the shoot apex were harvested from each plant and stored in RNAlater TissueProtect Tubes until transfer to a freezer at minus 80 degrees Celsius. RNA was isolated using a CTAB-based extraction protocol, and only RNA with an integrity number greater than 8.0 was used for library preparation. RNA from all 20 individuals was used for TagSeq, while RNA from one representative individual per species was additionally used for Iso-Seq library preparation. The authors note that because field-collected and potted plants were not phenologically synchronized, developmental stage cannot be regarded as controlled across species, and differences in expression should be interpreted primarily as species-associated patterns in this dataset, acknowledging that environmental and source effects cannot be fully disentangled.

Iso-Seq libraries were used solely to construct the reference transcriptome, with expression quantified separately from TagSeq libraries. For each species, one microgram of RNA from a single representative individual was used for cDNA synthesis and SMRTbell library construction. Species-specific barcoded double-stranded cDNAs were pooled at equal amounts and size-selected using a BluePippin, collecting 1 to 2 kilobase pair, 2 to 3 kilobase pair, 3 to 6 kilobase pair, and 5 to 10 kilobase pair fractions. The pooled barcoded library was sequenced on a single SMRT cell on a PacBio Sequel II system, yielding 3,599,918 circular consensus sequence reads totaling 8.11 gigabase pairs, with a mean read length of 2253 base pairs and median read quality of Q40. Demultiplexing assigned 848,223 reads to A. aethusifolius, 596,408 to A. dioicus, 928,092 to P. caespitosum, 639,089 to S. lucida, and 334,053 to S. splendens. TagSeq libraries prepared from all 20 individuals were sequenced single-end at 93 base pairs on an Illumina HiSeq 4000, with all libraries prepared in a single batch so that library preparation and sequencing batch were not confounded with species.

The assembled transcriptomes showed substantial variation among species. A. aethusifolius yielded the highest number of full-length reads at 847,960 and produced 59,679 nonredundant transcripts, whereas S. splendens, despite having the fewest full-length reads at 333,890, still produced 31,518 transcripts. Median transcript lengths ranged from 1607 to 2315 base pairs across species, and S. lucida contained the longest observed transcript at 8754 base pairs. Transcriptome completeness was assessed with BUSCO v5.8.2 against the embryophyta_odb10 lineage of 1614 BUSCO groups. Four of the five assemblies were closely comparable: A. dioicus and P. caespitosum each reached 82.8 percent completeness, A. aethusifolius 81.8 percent, and S. lucida 81.7 percent. S. splendens was appreciably less complete at 71.4 percent, with 26.6 percent of BUSCOs missing compared with 14.3 to 17.0 percent for the other species. This difference tracks sequencing depth, as S. splendens received approximately one-third the Iso-Seq yield obtained for P. caespitosum. Analyses requiring ortholog recovery in S. splendens are therefore depth-limited, and reciprocal-best-hit ortholog identification recovered a mean of 4107 pairs for comparisons involving S. splendens compared with 6336 for those that do not, a reduction of approximately 35 percent.

Exploratory analyses of global expression patterns revealed that both hierarchical clustering and principal component analysis separated samples into three major groups corresponding to the genera Aruncus, Spiraea, and Petrophytum, consistent with taxonomy and previously reported phylogenetic relationships within Spiraeeae. In the PCA, the first two principal components explained 54.94 percent and 33.75 percent of the total variance, respectively, and biological replicates clustered tightly by species. Notably, S. lucida and S. splendens grouped together despite being sampled under different conditions, with one potted and the other field-collected, indicating that growth condition was not the principal driver of this structure. The authors nonetheless cannot separate species identity from sampling condition, so they describe the observed patterns of differential expression as species-associated, without attributing them to genetically based regulatory divergence.

Pairwise comparisons among the five species yielded ten contrasts, applying thresholds of a log2 fold change of at least 2 and a false discovery rate below 0.01. Transcript abundance differed extensively among the species. In the comparison between A. aethusifolius and P. caespitosum, 6836 transcripts were significantly upregulated and 4442 were downregulated. The comparison between A. aethusifolius and S. splendens identified 5350 upregulated and 4135 downregulated transcripts, while the comparison between S. lucida and P. caespitosum yielded 3811 upregulated and 2755 downregulated transcripts. As expected, comparisons between more closely related species, such as A. aethusifolius versus A. dioicus and S. lucida versus S. splendens, resulted in substantially fewer differentially expressed transcripts. Analysis of transcription factor families among these transcripts revealed that the bHLH and bZIP families showed the greatest variability in differential expression counts, particularly in contrasts involving A. aethusifolius versus S. splendens and S. lucida versus P. caespitosum, with moderate variability for WRKY and C2H2 families. The authors caution that because expression was measured against a reference assembled from species-specific Iso-Seq contigs, reads mapped more efficiently to transcripts contributed by their own species, at rates of 70.5 to 85.1 percent, than to those contributed by others, at 25.1 to 77.0 percent, so the magnitude of between-species expression differences is expected to be inflated relative to true biological divergence, with cross-genus comparisons more affected than within-genus ones.

Functional enrichment analyses revealed several taxon-specific patterns. KEGG pathway enrichment showed that pathways related to aminoacyl-tRNA biosynthesis, fatty acid metabolism, and peroxisome function were enriched exclusively in comparisons involving A. aethusifolius versus P. caespitosum and A. aethusifolius versus S. splendens, while ribosome and oxidative phosphorylation pathways were enriched only in comparisons involving A. dioicus against S. lucida, S. splendens, and P. caespitosum. Gene Ontology enrichment for molecular function was dominated by binding and catalytic activity terms, including ion binding, metal ion binding, and ATP binding, with certain functions such as RNA binding showing comparison-specific enrichment. Cellular component analysis highlighted widespread involvement of intracellular organelles, with mitochondrial terms uniquely enriched in the S. lucida versus P. caespitosum comparison and several chloroplast-related terms enriched only in comparisons involving S. lucida against A. aethusifolius and S. splendens. The authors interpret these recurring organelle- and energy-associated categories as suggesting that lineage differences in cellular energy balance and stress-associated metabolic processes are a recurring component of transcriptomic differentiation among these taxa, while emphasizing that these enrichments do not by themselves establish adaptive or environmentally induced mechanisms.

Sequence-level analyses provided the clearest evolutionary signal in the study. The researchers analyzed distributions of synonymous substitution rates across orthologous gene pairs, which resolved into three tiers corresponding to the subclade structure of the tribe. The two congeneric comparisons showed the lowest divergence, with A. aethusifolius versus A. dioicus exhibiting the sharpest peak at a median Ks of 0.032 and S. splendens versus S. lucida next at a median Ks of 0.066. Intermediate distributions were recovered for the two intergeneric comparisons within the Spiraea–Petrophytum subclade, and the six comparisons between Aruncus and the remaining genera, which span the two subclades, showed the greatest divergence with median Ks values of 0.226 to 0.260. Across all ten comparisons, many orthologous gene pairs exhibited Ka/Ks ratios well below 1, indicating pervasive purifying selection on protein-coding sequences, which the authors describe as the dominant and most robust signal in these data. A small fraction of retained pairs showed Ka/Ks ratios greater than 1, ranging from 2.70 percent for A. dioicus versus P. caespitosum to 4.97 percent for A. aethusifolius versus A. dioicus, but these elevated-ratio pairs were concentrated at low synonymous divergence where estimates have high variance and are prone to inflation. The authors therefore interpret these loci as screening-level candidates for rapid sequence divergence, noting that confirming any individual candidate would require formal codon-based selection tests. Annotated candidates included heat shock proteins, metal ion transporters, zinc finger proteins, cytochrome P450 family members, and other diverse functional categories.

The study also screened the combined reference transcriptome for core leaf-developmental regulators, recovering members of the principal families that pattern leaf shape, size, and polarity, including CIN-type TCP transcription factors, GROWTH-REGULATING FACTORs with their GIF/AN3 coactivators, YABBY and HD-ZIP III/REVOLUTA genes, and BELL and WOX homeobox genes. The developmental-family orthologs that could be reliably tested and confidently annotated were uniformly under purifying selection, with Ka/Ks ratios of approximately 0.10 to 0.25, consistent with the genome-wide pattern and with the strong functional constraint expected for these conserved regulators. The authors state that the conservation of their protein-coding sequences indicates that the marked differences in leaf form among these species are unlikely to arise from amino acid changes in the regulators themselves, pointing instead toward regulatory or expression-level evolution as the more probable basis. Notably, the class I KNOX program, which promotes leaflet initiation in compound-leaved species, was essentially unrepresented, with only a single transcript recovered, so the compound leaves of Aruncus are not accompanied by a detectable class I KNOX signature in these data. Because leaf morphology is confounded with taxon in this sampling design and no quantitative morphometric measurements were made, the authors do not attempt to associate expression differences with leaf form, instead presenting the work as a comparative resource and a bounded set of candidate genes, pathways, and rapidly evolving loci to be tested with controlled developmental and ecological sampling, ideally in a common garden design.

Subject of Research: Comparative transcriptomics of five species in the Rosaceae tribe Spiraeeae

Article Title: Long‐Read Reference Transcriptomes for Five Species of Rosaceae Tribe Spiraeeae: A Resource for Investigating Sequence Divergence, Species‐Associated Expression Variation, and Leaf‐Developmental Candidate Genes

Article References: Mafakheri, A. M., & Potter, D. (2026). Long‐Read Reference Transcriptomes for Five Species of Rosaceae Tribe Spiraeeae: A Resource for Investigating Sequence Divergence, Species‐Associated Expression Variation, and Leaf‐Developmental Candidate Genes. Plant Direct, 10(11), Article e70202. https://doi.org/10.1002/pld3.70202

Image Credits: AI Generated

DOI: 10.1002/pld3.70202

Keywords: transcriptomics, Rosaceae, Spiraeeae, Iso-Seq, differential expression, Ka/Ks, purifying selection, leaf development, nonmodel species, BUSCO, WGCNA, orthologs